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A  model  is  developed  for  thermomechanical  behavior  of  defective,  low-symmetry  ceramic 
crystals  such  as  a-corundum.  Kinematics  resolved  are  nonlinear  elastic  deformation,  ther¬ 
mal  expansion,  dislocation  glide,  mechanical  twinning,  and  residual  lattice  strains  associ¬ 
ated  with  eigenstress  fields  of  defects  such  as  dislocations  and  stacking  faults.  Multiscale 
concepts  are  applied  to  describe  effects  of  twinning  on  effective  thermoelastic  properties. 
Glide  and  twinning  are  thermodynamically  irreversible,  while  free  energy  accumulates 
with  geometrically  necessary  dislocations  associated  with  strain  and  rotation  gradients, 
statistically  stored  dislocations,  and  twin  boundaries.  The  model  is  applied  to  describe  sin¬ 
gle  crystals  of  corundum.  Hardening  behaviors  of  glide  and  twin  systems  from  the  total 
density  of  dislocations  accumulated  during  basal  slip  are  quantified  for  pure  and  doped 
corundum  crystals.  Residual  lattice  expansion  is  predicted  from  nonlinear  elasticity  and 
dislocation  line  and  stacking  fault  energies. 
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1.  Introduction 

By  standard  definition,  a  ceramic  crystal  is  an  inorganic,  non-metallic,  ordered  solid.  Electronic  structures  of  ceramics  fea¬ 
ture  ionic  and/or  covalent  bonds  rather  than  metallic  bonds  found  in  ductile  metals.  The  nature  of  inter-atomic  forces  in 
ceramics  often  correlates  with  a  high  Peierls  barrier  (Peierls,  1940;  Nabarro,  1947;  Friedel,  1964),  inhibiting  dislocation  mo¬ 
tion  at  low  temperatures  and  leading  to  brittleness,  i.e.,  a  tendency  towards  fracture  over  slip. 

Anisotropic  and  non-cubic  crystals  are  of  particular  interest  in  the  present  work.  In  certain  hexagonal  ceramics  such  as 
silicon  carbide  (Zhang  et  al.,  2005),  basal  slip  is  often  the  preferred  inelastic  deformation  mechanism,  with  slip  resistances 
extremely  high  In  directions  normal  to  the  basal  plane.  This  phenomenon  occurs  similarly  In  many  hexagonal  metals,  includ¬ 
ing  certain  alloys  of  zirconium  (Tome  et  al.,  1991a:  McCabe  et  al.,  2009),  titanium  (Schoenfeld  and  Kad,  2002;  Mayeur  and 
McDowell,  2007),  and  magnesium  (Staroselsky  and  Anand,  2003;  Neil  and  Agnew,  2009),  though  in  some  cases  pyramidal 
slip  modes  are  possible.  Deformation  twinning,  as  opposed  to  slip,  is  often  the  only  viable  mechanism  for  accommodating 
deformation  normal  to  the  basal  plane.  In  lieu  of  fracture.  Twinning  is  often  favored  over  slip  in  cubic  metals  with  low  stack¬ 
ing  fault  energies  (Christian  and  Mahajan,  1995;  Kalldindi,  1998,2001);  twinning  may  also  influence  shear  band  formation  In 
cubic  metals  (Paul  et  al.,  2009).  In  addition  to  silicon  carbide,  other  low-symmetry  ceramics  of  interest  for  high  rate  appli¬ 
cations,  because  of  their  hardness  and  dynamic  compressive  strength,  include  titanium  diboride  (Bourne  and  Gray,  2002) 
and  alumina  (Bourne,  2006;  Bourne  et  al.,  2007). 

At  low  temperatures  and  pressures,  the  stable  phase  of  single  crystalline  alumina  (AI2O3)  is  a-corundum.  Although  other 
phases  exist  (Holm  et  al.,  1999),  henceforth  the  term  corundum  will  refer  only  to  the  a  phase.  The  Bravals  lattice  is 
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rhombohedral  (i.e.,  trigonal),  though  the  unit  cell  is  often  described  via  hexagonal  notation  (Kronberg,  1957;  Snow  and  Heu- 
er,  1973).  Corundum  is  centrosymmetric  and  hence  not  piezoelectric.  Corundum  that  is  red  in  color,  for  example  resulting 
from  chromium  doping,  is  commonly  called  ruby  (Chang,  1960;  Klassen-Neklyudova  et  al.,  1970;  Inkson,  2000).  Corundum  of 
all  other  colors  Is  often  referred  to  as  sapphire,  with  blue  sapphire  containing  trace  amounts  of  cobalt,  titanium,  or  iron,  for 
example.  Corundum  Is  extremely  hard,  with  a  value  of  9  on  Mohr’s  scale,  and  exhibits  a  very  high  Hugoniot  Elastic  Limit 
(HEL,  the  usual  measure  of  dynamic  yield  strength  under  uniaxial  strain  conditions,  at  high  pressure),  with  values  In  excess 
of  20  GPa  for  single  crystals  of  certain  orientations  (Graham  and  Brooks,  1971).  Single  crystals  of  adequate  purity  are  trans¬ 
parent,  with  uses  In  optics  and  electronics.  However,  the  HEL  and  toughness  of  commercial  grade,  polycrystalllne  alumina 
are  lower  than  that  of  single  crystalline  corundum,  often  thought  a  result  of  porosity  and  amorphous  phases  in  the  vicinity  of 
grain  boundaries  (Bourne  et  al.,  2007).  Industrial-grade  polycrystalline  alumina  manufactured  by  sintering  is  generally  not 
transparent. 

The  present  study  focuses  on  mechanisms  of  elasticity,  plastic  slip,  and  deformation  twinning.  Each  mechanism  is  ad¬ 
dressed  independently  via  an  individual  term  within  a  three-term  multiplicative  decomposition  of  the  deformation  gradient 
(Kratochvil,  1972;  Clayton  et  al.,  2005).  According  to  this  representation,  plastic  deformation  is  deemed  lattice-preserving 
(i.e.,  dislocation  glide  does  not  affect  the  lattice  vectors  or  stored  elastic  energy  of  crystal),  following  Bilby  et  al.  (1957), 
who  introduced  a  two-term  multiplicative  decomposition  for  describing  elasticity  and  plasticity  of  crystalline  solids  sub¬ 
jected  to  large  deformations.  In  contrast,  twinning  is  modeled  here  distinctly  from  dislocation  plasticity  via  the  use  of  an  iso- 
choric  intermediate  term  in  the  deformation  gradient  decomposition,  following  Kratochvil  (1972),  Clayton  et  al.  (2005),  who 
used  Intermediate  terms  within  three-term  decompositions  to  account  for  irreversible  deformations  that  are  not  lattice-pre- 
serving.  Microscopic  strain  fields  associated  with  defects  in  the  lattice  may  also  contribute  to  the  intermediate  deformation 
mapping.  Here,  isotropic  contributions  of  the  local  fields  follow  from  multiscale  volume  averaging  and  nonlinear  elastic  anal¬ 
ysis  of  self-equilibrated  bodies  with  defects  (Seeger  and  Haasen,  1958;  Toupin  and  Rivlin,  1960;  Teodosiu,  1982;  Wright, 
1982;  Clayton  and  Bammann,  2009).  following  the  theory  of  continuum  crystal  plasticity  (Hutchinson,  1976;  Teodosiu 
and  Sidoroff,  1976;  Asaro,  1983),  plastic  deformation  takes  place  via  slip  on  one  or  more  discrete  systems.  Twinning  takes 
place  via  energy  Invariant  shears  of  predefined  magnitude,  with  the  rate  of  shearing  determined  by  the  rate  of  increase  in 
volume  fraction  of  the  twin  relative  to  the  parent  (Chin  et  al.,  1969;  Van  Houtte,  1978;  Staroselsky  and  Anand,  2003).  Ther¬ 
modynamically  reversible  deformation  of  the  lattice  is  addressed  via  the  elastic  term  in  the  decomposition,  encompassing 
recoverable  deformation  associated  with  mechanical  stress  and  stress-free  thermal  expansion/contraction.  Multiscale  aver¬ 
aging  concepts  are  invoked  to  describe  effects  of  twinning  on  the  thermoelastic  response,  leading  to  effective  anisotropic 
elastic  coefficients  and  dissipation  rates  that  evolve  In  conjunction  with  the  volume  fractions  of  twins  in  an  element  of  fixed 
mass  in  the  crystal. 

A  constitutive  framework  based  on  Internal  state  variable  theory  provides  thermodynamic  relationships  among  indepen¬ 
dent  and  dependent  state  variables  as  well  as  appropriate  driving  forces  for  evolution  of  internal  variables  and  rates  of 
Inelastic  deformations.  Here,  the  former  consist  of  dislocation  and  twin  boundary  densities,  and  the  latter  Include  rates  of 
slip  and  of  twinned  volumes.  Dislocation  densities  consist  of  geometrically  necessary  dislocations  associated  with  slip  gra¬ 
dients  (Nye,  1953;  Ashby,  1970;  fleck  et  al.,  1994;  Arsenlis  and  Parks,  1999;  Rezvanian  et  al.,  2007)  and  statistically  stored 
dislocations  associated  with  homogeneous  plastic  flow  and  dislocation  loops  (Ashby,  1970;  Arsenlis  and  Parks,  1999;  Rez¬ 
vanian  et  al.,  2007).  Interface  partial  dislocations  at  propagating  twin  boundaries  (Scott  and  Orr,  1983)  are  demonstrated 
to  also  contribute  to  the  geometrically  necessary  dislocation  density  tensor,  following  kinematic  arguments  regarding  inelas¬ 
tic  deformation  incompatibility  in  plastically  deformed  single  crystals  (Clayton  et  al.,  2004a,b).  Thermodynamic  restrictions 
on  kinetic  relations  follow  naturally  from  energy  conservation  requirements  and  the  entropy  inequality  (Eckart,  1948;  Cole¬ 
man  and  Curtin,  1967;  Teodosiu,  1970;  McDowell,  2005).  Defects  In  corundum  are  known  to  contribute  to  lattice  curvature 
(Nye,  1953)  and  work  hardening  (Klassen-Neklyudova  et  al.,  1970;  Pletka  et  al.,  1977,  1982). 

The  theoretical  model  formalized  here  can  be  used  to  represent  corundum  single  crystals  within  an  aggregate  of  polycrys¬ 
talllne  alumina,  for  example  In  mesoscale  numerical  simulations  (Bourne,  2006)  of  grain  interactions,  stress  wave  propaga¬ 
tion,  and  defects.  Additionally,  the  model  may  enable  refinement  of  macroscopic  continuum  models  used  in  engineering 
design  (Rajendran,  1994),  providing  insight  into  important  mechanisms  occurring  at  the  scale  of  Individual  grains. 

The  present  application  focuses  on  alumina  because  of  the  availability  of  experimental  data  quantifying  the  mechanical 
(I.e.,  stress-deformation)  response  of  single  crystals  of  various  orientations  across  a  wide  range  of  temperatures  and  strain 
rates  (Graham  and  Brooks,  1971 ;  Tymiak  and  Gerberich,  2007;  Rodriguez  et  al.,  2008).  While  various  kinetic  relationships  for 
slip  (Lagerlof  et  al.,  1994)  and  twinning  (Scott  and  Orr,  1983)  have  been  posited  to  describe  the  results  of  specific  experi¬ 
ments,  a  more  general  continuum  model  Is  needed  to  collectively  explain  the  behavior  of  corundum  over  a  range  of  loading 
conditions. 

In  a  previous  paper  (Clayton,  2009),  a  theoretical  model  was  developed  to  address  anisotropic  nonlinear  elasticity,  slip, 
and  twinning  in  pure  sapphire.  The  focus  of  that  investigation  centered  on  the  shock  response  of  oriented  single  crystals  sub¬ 
jected  to  uniaxial  strain,  with  results  of  nonlinear  elastic  calculations  combined  with  experimental  HEL  data  to  provide 
bounds  on  critical  stresses  for  dislocation  glide  and  twin  nucleation.  Here,  that  theoretical  model  Is  extended  to  delineate 
geometrically  necessary  dislocations  resulting  from  gradients  of  inelastic  stretch  and  rotation  (i.e.,  deformation  incompati¬ 
bility)  as  well  as  statistically  stored  dislocations  associated  with  dislocation  dipoles  and  loops,  for  example.  Inclusion  of  geo¬ 
metrically  necessary  dislocations  requires  introduction  of  a  length  scale  (Eleck  et  al.,  1994;  Regueiro  et  al.,  2002)  In  the 
thermodynamic  framework  and  renders  the  model  nonlocal  according  to  many  labeling  schemes.  The  application  In  the 
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present  paper  focuses  on  strain  hardening,  defect  accumulation,  and  stored  energy  of  defects  in  doped  corundum,  in  addition 
to  the  pure  sapphire  considered  previously.  A  more  extensive  treatment  of  volume  changes  resulting  from  eigenstress  fields 
of  defects  is  also  given  here,  including  predictions  for  volume  changes  resulting  from  twin  boundary  energies. 

This  paper  is  organized  as  follows.  In  Section  2,  physical  descriptions  of  elasticity,  plasticity,  and  twinning  are  provided,  in 
a  more  extensive  treatment  than  given  previously  (Clayton,  2009).  The  descriptions  serve  to  distinguish  among  the  three 
deformation  mechanisms  and  provide  sufficient  physical  basis  for  the  corresponding  theoretical  framework  that  follows 
in  Section  3.  The  framework  of  Section  3  is  generic  enough  to  apply  to  any  ciystalline  solid  that  undergoes  large  deformations 
via  elasticity,  plasticity,  and/or  twinning,  including  metallic  and  ceramic  crystals.  In  Section  4,  the  model  is  specialized  to 
describe  the  behavior  of  corundum  single  crystals.  Section  4  begins  with  a  review  of  the  crystal  structure  of  corundum 
(i.e.,  atomic  positions,  lattice  parameters,  and  slip  and  twinning  systems)  and  its  fundamental  physical  properties.  Hardening 
behaviors  of  glide  and  twin  systems  from  the  density  of  dislocations  accumulated  during  basal  slip  observed  in  experiments 
on  pure  alumina  (Pletka  et  al.,  1977;  Castaing  et  al.,  2002)  and  corundum  doped  with  Cr,  Mg,  orTi  cations  (Pletka  et  al.,  1982) 
are  addressed.  Stored  residual  energies  associated  with  stress  fields  of  dislocations,  stacking  faults,  and  twin  boundaries  are 
quantified.  Residual  elastic  volume  changes  in  pure  and  doped  corundum  are  predicted  from  nonlinear  elasticity  theory  and 
dislocation  line  and  stacking  fault  energies.  New  comparisons  of  model  predictions  for  rhombohedral  twinning  with  exper¬ 
iments  (Scott  and  Orr,  1983;  Castaing  et  al.,  2002)  are  provided.  Finally,  comparisons  of  model  features  (e.g.,  basal  glide  resis¬ 
tance,  basal  dislocation  energy,  and  expansion  associated  with  basal  stacking  faults)  with  results  of  atomic  simulations 
performed  elsewhere  (Marinopoulos  and  Elsasser,  2001;  Bodur  et  al.,  2005;  Zhang  et  al.,  2007,  2008;  Nishimura  et  al., 
2009)  are  reported.  The  present  paper  is  self-contained;  hence,  some  overlap  of  material  in  Sections  3  and  4  with  content 
of  a  previous  paper  (Clayton,  2009)  is  inevitable  because  all  definitions  of  mathematical  symbols  and  values  of  material  prop¬ 
erties  entering  the  present,  more  refined  theoretical  model  and  application  are  included  here  for  completeness.  Differences 
in  theory  and  results  from  the  previous  paper  (Clayton,  2009)  are  highlighted  as  they  appear. 

The  following  notation  is  used.  Scalars  and  individual  components  of  vectors  and  tensors  are  written  in  italic  font,  while 
vectors  and  tensors  are  written  in  bold  font.  Einstein’s  summation  convention  applies  for  repeated  indices.  The  ■  symbol  de¬ 
notes  the  scalar  product  of  vectors  (a  ■  b  =  a“ba  =  a'b]  +  a^b2  +  a^bi),  while  ®  indicates  the  outer  product  ((a  ®  b)”*”  =  a°b*'). 
Juxtaposition  of  second-rank  tensors  implies  summation  over  one  set  of  adjacent  indices  ((AB)“  =  A“‘’6(,c)-  Summation  over 
two  sets  of  adjacent  indices  is  denoted;  (A  :  B  =  A“‘’Baj,).  Indices  in  parentheses  are  symmetric  (2A(ab)  =  +Aba)\  indices  in 
braces  are  skew  (2A[aj,]  =  Aab  -  Aba)-  Superposed  ■,  T,  and  -1  denote  material  time  differentiation,  transposition,  and  inver¬ 
sion,  respectively,  and  subscripted  commas  denote  partial  differentiation. 

2.  Background:  deformation  mechanisms 

2.1.  Elasticity 

A  crystal  is  said  to  deform  elastically  in  the  absence  of  generation  or  motion  of  defects.  At  the  atomic  scale,  elastic  defor¬ 
mation  alters  relative  distances  and/or  orientations  among  neighboring  atoms  within  each  crystallographic  unit  cell.  Result¬ 
ing  changes  in  inter-atomic  forces  produce  mechanical  stress  when  the  crystal  is  viewed  as  a  continuous  solid.  Removal  of 
mechanical  stresses  restores  the  original  inter-atomic  bond  lengths  and  angles  without  dissipation  of  energy;  hence,  elastic 
deformation  is  said  to  be  thermodynamically  reversible.  Here,  elastic  deformation  also  includes  changes  in  average  inter¬ 
atomic  bond  vectors  induced  by  changes  in  temperature.  Increases  in  thermal  energy,  i.e.,  local  atomic  vibrations,  usually 
correlate  with  expansion  of  the  lattice  in  the  absence  of  mechanical  stresses.  Such  thermal  deformations  are  deemed  revers¬ 
ible,  since  atomic  nuclei  will  return  to  their  mean  reference  positions  as  the  temperature  is  restored  to  its  original  value.  Let 

denote  the  two-point  tensor  of  possibly  large  elastic  deformation.  Let  ao  denote  the  vector  between  sites  in  a  conventional 
unit  cell  such  that  Uq  =  (ao  ■  ao)’^^  is  a  conventional  lattice  parameter.  Primitive  Bravais  lattice  vectors  and  conventional  lat¬ 
tice  vectors  deform  via  the  elastic  deformation  according  to  the  Cauchy-Born  hypothesis  (Born  and  Huang,  1954;  Ericksen, 
1984).  This  is  written  a  =  F^ao,  with  a  the  conventional  lattice  vector  in  the  deformed  ciystal.  An  elastic  deformation  by  sim¬ 
ple  shear  of  a  non-cubic  lattice  is  shown  in  Fig.  1(a).  In  the  present  description,  electric  polarization  and  relative  internal 


Fig.  1.  Deformation  mechanisms  and  lattice  vectors:  (a)  elasticity,  (b)  slip,  and  (c)  twinning. 
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shifts  between  sub-lattices  in  non-centrosymmetric  crystals  (Born  and  Huang,  1954;  Cousins,  1978)  are  not  addressed 
explicitly;  rather,  elastic  deformation  of  a  volume  element  is  uniform  over  all  inter-atomic  bond  vectors  within  that  element. 

2.2.  Plasticity 

Plastic  deformation  as  defined  here  takes  place  via  glide  of  dislocations  of  edge,  screw,  and/or  mixed  character,  including 
loops,  and  encompassing  cross-slip  but  not  climb,  the  latter  which  requires  further  generalization  (McDowell  and  Mooseb- 
rugger,  1992).  As  full  dislocations  travel  through  a  region  of  the  ciystal,  the  shape  of  the  material  will  change,  but  inter¬ 
atomic  distances  remain  the  same,  so  long  as  no  defects  are  left  behind  in  that  region.  In  this  sense,  plastic  deformation 
is  said  to  be  lattice-preserving  or  “lattice  invariant”  (Bilby  et  al.,  1957).  Mechanical  stresses  are  conventionally  required 
to  enact  the  net  glide  of  dislocations  (apart  from  random  thermal  fluctuations).  For  example,  resolved  shear  stresses  must 
exceed  the  Peierls  barrier  in  the  context  of  lattice  statics  or  Schmid’s  limit  in  the  context  of  continuum  slip  (Hirth  and  Lothe, 
1982).  Plastic  deformation  is  thermomechanically  irreversible,  since  the  reference  shape  of  the  material  is  not  recovered 
upon  removal  of  mechanical  stresses,  and  since  heat  is  dissipated  by  moving  dislocations  as  a  result  of  lattice  friction,  pho¬ 
non  drag,  and  other  mechanisms  (Kocks  et  al.,  1975;  Gilman,  1979).  Because  the  lattice  remains  unchanged  apart  from  steps 
on  the  surface  of  the  crystal,  plastic  deformation  itself  does  not  affect  the  strain  energy  of  the  crystal.  However,  defects  gen¬ 
erated  during  plastic  deformation  that  remain  within  the  material  lead  to  energy  storage  as  a  result  of  local  stress  fields  in¬ 
duced  by  these  defects  and  their  core  energies.  Let  F’’  denote  the  two-point  tensor  of  possibly  large  plastic  deformation. 
Plastic  deformation  by  simple  shear  is  illustrated  for  a  non-cubic  lattice  in  Fig.  1(b).  Lattice  vector  ao  is  unchanged  by  the 
plastic  deformation. 

2.3.  Twinning 

Deformation  twinning  results  in  two  connected  regions  in  the  lattice  separated  by  a  twin  boundary  (i.e.,  the  habit  plane  or 
dotted  line  on  the  right  side  of  Fig.  1(c))  whose  shape  deformations  differ  by  a  simple  shear.  The  original  lattice  is  termed  the 
“parent”,  while  the  sheared  lattice  is  termed  the  “twin”.  Atomic  positions,  and  hence  corresponding  bond  vectors  between 
these  atoms,  within  each  region  differ  by  a  finite  rotation,  typically  either  a  reflection  or  180“  rotation  (Christian  and  Maha- 
jan,  1995),  though  more  general  relationships  are  possible.  The  stacking  sequence  of  atomic  planes  in  the  twin  is  altered  with 
respect  to  that  in  the  parent,  and  hence  crystals  with  low  stacking  fault  energies  are  often  more  prone  to  twinning  than  those 
with  high  stacking  fault  energies.  Nucleation  and  propagation  of  deformation  twins  are  thought  to  take  place  by  one  or  more 
mechanisms,  often  involving  formation  and  motion  of  partial  dislocations  (e.g.,  dissociation  of  full  dislocations  into  partials) 
and  atomic  shuffles  sometimes  needed  to  maintain  orientation  relationships  between  twin  and  parent  (Bilby  and  Crocker, 
1965;  Zanzotto,  1996).  Thus,  twinning  is  not  regarded  as  lattice-preseiT/ing  in  the  sense  of  slip,  since  twinning  involves  rota¬ 
tion  of  the  lattice,  as  evidenced  by  texture  measurements  (Van  Houtte,  1978;  Tome  et  al.,  1991a).  Involvement  of  partial  dis¬ 
locations  and  atomic  shuffles  induces  the  rotational  change  of  the  twinned  lattice  relative  to  that  of  the  parent.  Deformation 
twinning  is  also  distinguished  from  plastic  slip  in  that  the  former  occurs  by  collective  motion  of  defects,  resulting  in  a  quan¬ 
tized  amount  of  shear  that  preseiT/es  the  particular  orientation  relationship  between  the  twin  and  parent.  In  contrast,  plastic 
deformation  may  result  in  shearing  of  any  magnitude,  with  the  lower  limit  of  relative  displacement  of  atomic  planes  asso¬ 
ciated  with  the  Burgers  vector  for  slip.  The  shear  strain  associated  with  twinning  is  deemed  mechanically  irreversible,  since 
twins  considered  here  remain  in  single  crystals  after  mechanical  stresses  are  removed.  Another  difference  between  slip  and 
twinning  is  that  twinning  is  polar  (i.e.,  unidirectional)  while  often  slip  is  not.  Lattice  geometry  precludes  twinning  shears  of 
equal  magnitude  and  opposite  directions  on  the  same  plane,  while  typically  slip  may  occur  in  opposite  directions  on  the 
same  plane,  though  resistances  to  slip  in  opposite  directions  on  the  same  plane  may  differ  (Lee  et  al.,  1999;  Xu  et  al.,  2004). 

Unstressed  twinned  regions  of  the  crystal  far  from  internal  boundaries  or  defects  possess  the  same  strain  energy  density 
as  the  unstressed  parent  Qames,  1981 ;  Zanzotto,  1996);  hence,  twinning  shears  are  said  to  be  energy  invariant.  However,  the 
energy  density  increases  relative  to  that  of  a  perfect  lattice  in  the  vicinity  of  twin  boundaries  (e.g.,  appropriate  stacking  fault 
energies).  The  mechanical  work  done  during  deformation  twinning  is  dissipative,  resulting  from  the  defect  motion  (e.g.,  par¬ 
tial  dislocation  glide)  associated  with  shearing.  Possible  energy  storage  is  associated  only  with  defects  left  behind  in  the  crys¬ 
tal,  for  example  those  comprising  the  twin  boundaiy.  From  continuum  thermomechanics  considerations,  the  driving  force 
for  twin  propagation  is  the  resolved  shear  stress  in  the  direction  of  twinning  shear.  The  resistance  to  deformation  twinning 
is  often  modeled  analogously  to  slip,  that  is,  twinning  proceeds  when  the  resolved  stress  attains  a  critical  value  that  may 
depend  on  temperature  (Lagerlof  et  al.,  1994;  Wu  et  al.,  2007).  With  accumulated  slip  and  deformation  twinning,  strain  hard¬ 
ening  of  the  crystal  may  take  place  via  interactions  among  different  twins  in  the  crystal,  interactions  among  different  slip 
systems,  and  interactions  between  mobile  dislocations  and  twins  (e.g.,  twin  boundaries  may  serve  as  barriers  to  dislocation 
glide).  Twins  may  also  nucleate  cracks  and  vice-versa  (Christian  and  Mahajan,  1995).  Detwinning,  i.e.,  restoration  of  the 
twinned  lattice  to  its  original  orientation,  is  physically  possible,  though  is  often  more  applicable  to  phase  transformation 
phenomena  (Bhattacharya,  1991 ;  Thamburaja  et  al.,  2009)  and  less  applicable  to  mechanical  twinning  in  the  context  of  mon¬ 
atomic  loading.  In  certain  metals  such  as  Mg  alloys,  however,  deformation-induced  detwinning  can  be  important,  especially 
during  load  sequences  involving  strain  path  changes  (Proust  et  al.,  2009). 

Twins  are  usually  classified  as  type  I,  type  II,  or  compound.  In  centrosymmetric  ciystals,  Bravais  lattice  vectors  in  the  twin 
and  parent  for  a  type  I  twin  are  related  by  either  a  reflection  in  the  habit  plane  or  rotation  of  180“  about  tbe  direction  normal 
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to  this  plane.  For  a  type  II  twin,  lattice  vectors  are  related  by  either  a  rotation  of  180°  about  the  shear  direction  or  a  reflection 
in  the  plane  normal  to  the  shear  direction.  In  crystals  with  a  center  of  symmetry,  the  rotation  (Van  Houtte,  1978;  Christian 
and  Mahajan,  1995) 

f  2mo(g)mo-l  (type  I), 

\2so®So-l  (type  II), 

relates  a  lattice  vector  in  the  parent,  ao,  to  a  vector  in  the  twin,  a,  via  a  =  Q,ao,  as  shown  on  the  right  of  Fig.  1(c).  The  unit 
normal  to  the  habit  plane  is  mo.  and  the  direction  of  shear  is  Sq.  The  second-order  unit  tensor  is  denoted  by  1.  For  centro- 
symmetric  crystals,  reflection  -Q,  is  often  used  instead  of  (1 )  to  describe  the  orientation  of  the  twin  relative  to  the  parent.  In 
crystals  lacking  a  center  of  symmetry,  the  two  rotations  listed  in  (1 )  are  augmented  by  their  negatives,  i.e.,  two  complemen¬ 
tary  reflection  operations  (Christian  and  Mahajan,  1995),  with  all  four  operations  then  crystallographically  distinct. 

3.  Continuum  theory:  nonlinear  elasticity,  slip,  and  twinning 

A  constitutive  framework  for  crystals  undergoing  large  thermoelastic,  plastic,  and  twinning  deformations  is  developed. 
The  framework  conforms  to  many  established  principles  of  continuum  mechanics  and  thermodynamics  of  single  crystal 
behavior,  as  described  in  Teodosiu  (1970),  Teodosiu  and  Sidoroff  (1976),  Clayton  (2005),  though  these  prior  works  did  not 
explicitly  consider  twinning  or  possible  residual  volume  changes  resulting  from  dislocations  and  twin  boundaries  or  stacking 
faults. 


3.1.  Kinematics 


A  small  volume  element  of  crystalline  material  of  fixed  mass  is  assigned  reference  coordinates  X.  Let  x  =  ip{X,  t)  denote 
spatial  coordinates  of  the  element,  with  tp  the  motion.  Deformation  gradient  F  for  the  element  is 


dx 

ax’ 


(2) 


decomposed  multiplicatively  into  a  series  of  terms: 

F  =  F^J’/^F'F'’  =  F’^Pf’’  =  F^F  =  F'-P.  (3) 

Here,  F®  accounts  for  recoverable  thermoelastic  deformation  and  rigid  body  rotation,  P  =  JV^p'  accounts  for  defect  kine¬ 
matics  that  alter  the  lattice,  and  F**  accounts  for  lattice-preserving  plastic  slip.  Twinning  is  modeled  by  term  F‘,  accounting 
for  the  average,  irreversible  shape  deformation  resulting  from  one  or  more  twins  that  may  nucleate  and  propagate  in  a  vol¬ 
ume  element  of  fixed  mass  of  crystal.  Volume  changes  associated  with  defects  are  addressed  by  scalar/.  The  total  irreversible 
deformation  is  F  =  F'F*”.  The  total  elastic  lattice  deformation  (recoverable  and  residual)  is  F*'  =  F^/V^l,  and  the  remaining 
deformation  from  defect  motion  (dislocation  glide  and  twinning)  is  F  =  F'F**. 

Denote  by  g(x)  and  G(X)  metric  tensors  associated  with  possibly  curvilinear  coordinate  systems  in  configurations  B  and 
Bo,  with  components  =  ga  ■  Sb  ^nd  Gab  =  Ga  ■  Cb  (Eringen,  1962).  Basis  vectors  in  spatial  and  referential  coordinate  sys¬ 
tems  are  g„  and  Ga,  respectively  (a,  A  =  1,2,3).  Let  all  intermediate  configurations  be  referred  to  an  external  Cartesian  coor¬ 
dinate  system  with  metric  tensor  components  i5a;)(a  =1,2, 3).  Let  the  Jacobian  determinants  of  deformation  mappings  in  (3) 
associated  with  volume  changes  be  formally  defined  as 

/=V^detF=/W=///,  /  =  \/d^detF^  /  =  detP,  /'’ =  y^detF'’.  (4) 

As  will  be  demonstrated  explicitly  later  in  (12),  the  isochoric  character  of  slip  and  twinning  leads  to  the  conditions 
/=/'’=!  (no  volume  changes  associated  with  twinning  or  slip)  and ]  =f]  =  /  (volume  changes  associated  only  with  recov¬ 
erable  and  residual  thermoelasticity).  The  residual  elastic  volume  change  /  arising  from  distributed  defects  such  as  disloca¬ 
tion  lines  is  derived  from  the  assertion  that  in  continuum  nonlinear  elasticity,  the  average  strain  of  a  body  containing 
residual  stress  fields  arising  from  internal  displacement  discontinuities  need  not  vanish  even  if  the  traction  on  its  external 
surfaces  vanishes  (Toupin  and  Rivlin,  I960:  Teodosiu,  1982;  Clayton  and  Bammann,  2009).  Presently,  only  volume  changes 
are  considered,  e.g.,  corresponding  to  random  defect  distributions  imparting  no  preferred  directions  in  average  residual  elas¬ 
tic  strains,  though  more  general  treatments  allowing  for  shape  changes  resulting  from  residual  stresses  associated  with  crys¬ 
tal  defects  at  multiple  length  scales  have  been  suggested  (Clayton  and  McDowell,  2003;  Clayton  et  al.,  2004a,  2005,  2006; 
McDowell,  2008;  Clayton  and  Bammann,  2009).  In  Section  4.4,  relationships  between/,  the  line  density  of  dislocations, 
and  nonlinear  thermoelastic  properties  are  given,  with  explicit  formulae  for/  listed  in  (65)  and  (66)  following  from  previous 
studies  (Zener,  1942;  Seeger  and  Haasen,  1958;  Toupin  and  Rivlin,  I960:  Holder  and  Granato,  1969;  Teodosiu,  1982;  Clayton 
and  Bammann,  2009).  In  the  present  context,  F**  and  F'  account  for  isochoric  deformations  resulting  from  respective  motion 
of  slip  dislocations  and  twinning  partials,  while/  accounts  for  volume  changes  resulting  from  the  residual  stress  fields  of  the 
defects  themselves. 
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Introduced  next  are  sets  of  contravariant  and  covariant  vectors  denoting  directions  and  planes,  respectively,  for  slip  and 
twinning.  When  referred  to  the  reference  lattice  prior  to  any  reorientation  by  twinning,  these  are  denoted  by  for 

each  slip  system  i,  and  |sj,,inj,|  for  each  twin  system  j.  The  total  number  of  slip  systems  is  n,  and  the  total  number  of  twin 
systems  is  w.  Reference  shearing  directions  and  plane  normals  are  all  of  unit  length,  and  each  pair  of  contravariant  shear 
direction  and  covariant  plane  normal  is  orthogonal: 

si,-m|,  =  0,  |s|,|  =  |m',|  =  1  (Vi=  s{,-mi,  =  0,  |si,|  =  |mi,|  =  l(y;  =  1, . . .  w).  (5) 

During  the  course  of  twinning,  one  or  more  parts  (i.e.,  twins)  of  the  volume  element  of  crystal  undergoes  a  rotation  rel¬ 
ative  to  the  parent.  In  a  volume  fraction  of  the  crystal  undergoing  twinning  via  mode  j,  slip  directions  and  slip  plane  normals 
transform  in  the  reference  configuration  according  to  the  usual  rules  for  contravariant  and  covariant  vectors,  that  is 

sjjj  =  Q/s|,,  m[,j  =  moQ'^  (Vi  =  1, ...  n  slip  systems ;Vj  =  1, . . .  w  twin  volumes),  (6) 

where  Qj  is  the  rotation  found  from  (1)  corresponding  to  particular  twin  system  j.  For  example,  if  j  is  a  type  1  twin, 
QJ  =  2mJ|  0  mj,  -  1,  while  ifj  is  a  type  II  twin,  Q/  =  2s|)  0  sj,  -  1.  Notice  from  (6)  that,  within  each  twinned  volume,  updated 
slip  directions  sj,j.  and  slip  plane  normals  mj|j  remain  orthogonal  and  of  unit  length  for  each  /.  For  simplicity,  successive  twin¬ 
ning  is  not  considered.  Flence,  secondary  twins  that  could  form  within  already  twinned  regions,  leading  to  reorientation  of 
the  twinning  systems  {sjj.mJ,},  are  not  represented.  Further  reorientation  of  transformed  directors  in  (6)  is  likewise  prohib¬ 
ited  once  the  twin  is  fully  formed.  Rotation  (6)  does  not  apply  to  the  volume  fraction  of  the  grain  comprising  the  parent. 
Plastic  deformation  F’’  and  twinning  deformation  F'  do  not  directly  alter  the  directions  associated  with  slip  and  twinning. 
The  former  is  lattice-preserving,  as  discussed  in  Section  2.2,  while  the  latter  affects  the  lattice  orientation  indirectly  via 
(6)  and  evolution  of  the  twin  volume  fraction  for  each  twin  system  j  to  be  discussed  later.  However,  thermoelastic  deforma¬ 
tion  and  residual  elastic  volume  changes  both  affect  the  lattice  directors: 


II 

o  — 

m'  =mi,F‘^ 

sj,  mj,  =  (Vi  =  1, . . 

.  n  e  parent) , 

s'  —  F^'s' 

=>or 

>  Soj  =J’^"sor  =  (Vi=l, 

.  ..n  e  twins  j  =  1 , . . . w) , 

(7) 

s'"  =  F^-sj), 

m'  ^m'gF*^ 

sj,  =  mj,  =  (Vj  =  1, . . 

.  w  e  parent) . 

The  effect  of  J  on  the  slip  directors  and  slip  plane  normals  was  omitted  in  earlier  work  (Clayton,  2009). 

An  additional  remark  on  notation  is  in  order.  Subscripts  and  superscripts  following  mathematical  objects  in  bold  font  cor¬ 
respond  to  slip  or  twin  systems  (e.g.,  i  or  j),  and  are  not  subject  to  usual  conventions  associated  with  the  index  notation  such 
as  Einstein’s  summation  convention.  Furthermore,  such  subscripts  and  superscripts  do  not  require  placeholder  periods  fre¬ 
quently  used  for  components  of  mixed  contra-covariant  tensors  (Schouten,  1954).  For  example,  the  deformation  gradient  is 
often  written  in  geometric  settings  as  F  =  F‘‘^ga  0  (Clayton  et  al.,  2005),  while  slip  direction  vector  for  system  i  within  twin 
system  j  is  written  here  as  sj  =  sj“g„  and  not  as  sj,  since  i  and  j  do  not  refer  to  components  of  a  vector  or  tensor  in  this  case. 

Fig.  2  depicts  the  physics  underlying  (2)-(7)  for  a  crystal  with  a  single  slip  system  and  a  single  twin.  Multiplicative 
decomposition  (3)  implies  a  series  of  configurations  of  the  material  element.  The  reference  configuration  is  labeled  Bo  with 


a 


b 


Fig.  2.  Deformations  and  slip  system  geometry  for  crystal  deforming  by  elasticity,  slip,  and  twinning:  (a)  total  deformation  and  (b)  residual  deformation. 
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corresponding  coordinates  X,  the  spatial  configuration  is  labeled  B  with  corresponding  coordinates  x,  and  the  elastically  un¬ 
loaded  intermediate  configuration  is  labeled  B.  Since  '  and  F  are  In  general  not  Integrable,  continuous  coordinates  span¬ 
ning  B  do  not  exist  (Clayton  et  al.,  2004b,  2005).  However,  elastic  and  Inelastic  deformations  act  as  tangent  maps  (Marsden 
and  Hughes,  1983)  between  configurations  via  :  TB  TB  and  F  :  TBg  TB.  Additional  configurations  Bp  and  B  are  im¬ 
plied  by  the  inelastic  tangent  maps  F** :  TBo  TBp.F'  :  TBp  TB,  andp/^1  :  TB  TB.  The  total  deformation  gradient  is 
illustrated  in  Fig.  2(a).  In  Fig.  2(b),  individual  maps  comprising  the  residual  deformation  F  =  F^^’F  =J’/3f'f'’  are  shown. 
For  clarity,  the  implicit  effect  of  twinning  on  the  orientation  of  the  plastically  slipped  lattice  in  configuration  Bp  is  not  shown; 
i.e.,  only  the  parent  Is  shown  in  the  illustration  of  configuration  Bp. 

The  spatial  velocity  gradient  following  from  (2)  and  (3)  is 


L  =  FF-’  =  F®F®-’  +  F'=F'F*-’F®-’  +  f'=F'F'’F‘’-’F'-’F'=-’  +  (1/3)))-’!. 

(8) 

The  Inelastic  velocity  gradient  referred  to  configuration  B  is 

FF-'  =F'f'-’  +F'F'’F'’-’f'-’  +(l/3)I/-’l  =L'  +  L'’  +  (l/3)i/-’l, 

(9) 

where 

L'  =  F'F'-’=  5]"/ys'o®m'o 

(10) 

twinning  in  parent  crystal 

results  from  twinning  shears,  and 


r 


=  f' 


slip  in  parent  crystal 


slip  in  twinned  crystal 


(11) 


results  from  slip  In  parent  and  twinned  domains.  In  (10),  ■'/  Is  the  predefined  shear  associated  with  twin  system),  a  positive 
scalar  that  is  fixed  for  all  twins  in  a  given  family  of  twin  systems,  in  (11 ),  y'  is  the  slip  rate  on  system  i  in  the  parent  grain,  and 
yj  is  the  slip  rate  on  system  i  within  reoriented  twin  fraction).  The  volume  fraction  of  ciystal  occupied  by  twin),  measured 
per  unit  volume  in  configuration  B,  is  labeled  by  the  scalar/1  ^  o,  with  time  rate /I  (Chin  et  al.,  1969;  Van  Houtte,  1978).  The 
total  volume  fraction  of  twinned  ciystal  is  fp  =  J2P’  subject  to  the  restriction  0  C  /r  ^  1-  Detwinning  is  not  considered; 
hence, /I  >  0.  In  the  Interior  summation  within  the  rightmost  term  of  (11),  the  slip  directors  and  slip  plane  normals  In 
the  twinned  regions  are  found  according  to  (6),  where  the  particular  form  of  Qj  corresponds  to  the  twin  with  associated  va¬ 
lue  of /I  In  the  outer  summation. 

Since  for  each  slip  or  twin  system,  the  shear  directions  and  plane  normals  are  orthogonal, 

trL'  =  trr’  =  0^]'=/  =  0^/ =)'’  =  !  {t  >  0),  (12) 

since  at  t  =  0,  F'  and  F*”  both  reduce  to  the  unit  tensor  in  corresponding  coordinate  systems  leading  to 
/{t  =  0)  =  f/t  =  0)  =  1,  and  where  trA  =  A“  for  a  second-order  matrix  A.  A  formal  derivation  of  (12)  Is  given  in  Appendix 
A.  Thus,  (9)-(ll)  properly  reflect  the  isochorlc  nature  of  slip  and  twinning,  and  together  with  (3)  and  (4),  require  that  all 
volume  changes  be  accommodated  thermoelastlcally  via  f  and/or  by  defect  generation  via),  such  that  the  first  of  (4)  reduces 
to)  =)")=)'. 

When )  =  1  and p  =  0  for  all),  F'  =  l.F  =  F*”,  and  L**  In  (11)  reduces  to  its  usual  definition  from  crystal  plasticity  theory 
(Teodosiu  and  Sidoroff,  1976;  Asaro,  1983;  Clayton,  2005).  In  this  simplified  case,  all  inelastic  deformation  occurs  via  slip,  F  is 
lattice-preserving,  and  the  lattice  directors  remain  unchanged  between  configurations  Bo  and  B.  On  the  other  hand,  when 
twinning  does  takes  place,  F'  does  not  act  as  a  true  elastic  lattice  deformation  in  the  sense  of  (7)  and  Born  and  Huang 
(1954),  since  only  that  part  of  the  lattice  within  the  twinned  volume  undergoes  a  transformation,  and  this  transformation 
occurs  via  rotation  (6)  and  does  not  include  any  stretch  of  the  lattice  directors. 

Remarks  on  the  order  of  terms  in  multiplicative  decomposition  (3)  are  now  merited.  The  recoverable  thermoelastic  defor¬ 
mation  F^  is  logically  placed  first  In  the  decomposition,  following  the  usual  scheme  of  crystal  plasticity  theory  (Teodosiu  and 
Sidoroff,  1976;  Asaro,  1983),  so  that  thermoelastic  unloading  proceeds  via  pre-multiplication  of  F  by  F®  ’  and  so  that  the 
lattice  director  vectors  are  mapped  to  the  current  configuration  via  (7).  The  remaining  three  terms  in  the  decomposition 
are  placed  in  order  of  decreasing  effect  on  the  crystal  lattice  (i.e.,  on  the  lattice  director  vectors).  TermJV^l  is  placed  second 
In  the  decomposition  because  it  affects  the  volume  of  the  crystal  and  imparts  a  stretch  to  the  lattice  directors  as  indicated  in 
(7);  however,  its  precise  placement  Is  of  apparently  lesser  mathematical  importance  because  this  term  is  spherical  (i.e.,  iso¬ 
tropic).  Term  F'  is  placed  next  because  twinning  affects  the  lattice  director  vectors  implicitly,  via  the  rotation  or  reflection 
operation  in  (6).  Term  F’’  is  placed  last  in  the  decomposition  because  plastic  slip  does  not  directly  affect  the  lattice  directors 
(Asaro,  1983).  Placement  of  F'  before  F**  in  (3)  does  not  imply  that  during  the  time  history  of  a  given  deformation  process, 
twinning  always  takes  place  before  slip,  or  vice-versa.  All  elastic  and  inelastic  deformation  modes  In  (3)  can  occur  simulta¬ 
neously  in  time.  This  argument  applies  similarly  in  the  context  of  traditional  crystal  plasticity  (i.e.,  no  twinning)  wherein 
F  =  F^F*”  holds:  for  example,  elastic  deformation  and  plastic  deformation  can  both  occur  simultaneously  at  a  given  time 
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as  an  increasing  load  is  applied  to  a  strain  hardening  crystal.  An  alternative  decomposition  in  which  F'  and  F’’  are  inter¬ 
changed  in  (3)  is  explored  in  Appendix  B. 

Defect  content  associated  with  incompatibility  of  the  total  lattice  deformation  F*'  (Teodosiu,  1970;  Regueiro  et  al.,  2002; 
Clayton  et  al.,  2004a, b,  2005)  is  measured  by  two-point  geometrically  necessary  dislocation  tensor  ac  satisfying 


B"  =  - 


f  =  -  /  F'^^dX''  =  f 

Jc  Jc  Ja 


afNAdA,  af  = 


[B.q’ 


(13) 


where  are  components  (indices  referred  to  configuration  B)  of  a  total  Burgers  vector  associated  with  circuit  c  in  the  spatial 
configuration  or  circuit  C  in  the  reference  configuration,  A  is  the  area  enclosed  by  C  with  unit  normal  components  Na,  and 
are  contravariant  components  of  the  permutation  tensor.  Stokes’s  theorem  is  used  to  convert  from  the  line  integral  to  the 
area  integral  in  the  third  equality  of  (13).  From  (3), 


_  p4BCpix  _  cABCplUpP/f  ,  ^BCpla  pP? 

“c  —  '  [B,q  ^  ‘  ‘  ll,C‘  B  ■ 


slip  gradients  twin  gradients 


(14) 


Contributions  of  plastic  slip  gradients  to  the  geometrically  necessaiy  dislocation  density  tensor  corresponding  to  the  first 
term  following  the  second  equality  in  (14)  are  well-documented  in  papers  dealing  with  the  continuum  theory  of  dislocations 
(Ashby,  1970;  Teodosiu,  1970;  Fleck  et  al.,  1994;  Arsenlis  and  Parks,  1999;  Voyiadjis  and  Abu  Al-Rub,  2007;  Rezvanian  et  al., 
2007;  Clayton  et  al.,  2008).  The  second  term  in  the  sum  on  the  right  arises  from  material  gradients  of  twinning  shear,  for 
example  gradients  of  twin  volume  fractions  arising  during  propagation  of  tapered  twins  (Scott  and  Orr,  1983);  the  contribu¬ 
tion  to  the  dislocation  density  tensor  in  this  case  would  be  partial  dislocations  at  interfaces  between  twin  and  parent  or  be¬ 
tween  intersecting  twins.  Applying  Nanson’s  formula  to  (13),  the  geometrically  necessary  dislocation  tensor  in  intermediate 
configuration  B  of  Fig.  2(b),  labeled  a,  is 


(15) 


where  J  =  ff  =  1  by  (12),  and  the  second  equality  in  (15)  reduces  to  Nye’s  relation  (Nye,  1953)  between  the  dislocation  den¬ 
sity  tensor  and  the  lattice  curvature  in  the  limit  of  small  elastic  deformations.  The  sum  on  the  far  right  of  ( 1 5 )  is  invoked  over 
dislocation  populations  of  index  label  i,  where  p'  is  the  length  per  unit  volume  of  dislocation  segments  with  unit  tangent  line 
vector  gi‘  and  Burgers  vector  b^.  Nye’s  original  treatment  (Nye,  1953)  was  in  part  based  on  observations  of  slip  traces  in  trans¬ 
parent  corundum,  and  hence  is  of  particular  relevance  in  the  present  application  to  ceramics  and  corundum  in  particular.  In 
(13)-(15),  sufficient  differentiability  of  F^.fVf’’,  andj  has  been  assumed.  Discontinuities  in  the  deformation  gradient  or  in 
the  Bravais  lattice,  for  example  singularities  across  slipped  regions  or  at  dislocation  cores  (Teodosiu,  1970, 1982)  and  jumps 
in  deformation  gradient  across  twin  boundaries  Qames,  1981 ;  Bhattacharya,  1991 )  are  not  resolved  explicitly  in  the  present 
continuum  framework  that  addresses  defects  via  continuous  distributions  (Willis,  1967). 

Non-dimensional  internal  state  variables  are  introduced  to  represent  energetic  changes  associated  with  two  other  kinds 
of  defects.  The  first  such  internal  state  variable  is  a  measure  of  the  density  of  statistically  stored  dislocations  (Ashby,  1970; 
Arsenlis  and  Parks,  1999;  Bammann,  2001;  Clayton  et  al.,  2004a,  2006)  that  accumulate  with  homogeneous  slip,  f  =  b^/p^, 
where  b  is  a  scalar  Burgers  vector— or  a  constant  on  the  order  of  a  lattice  parameter  when  the  crystal  exhibits  slip  on  systems 
with  different  Burgers  vectors— and  p^  is  the  total  length  of  such  dislocations  per  unit  volume  in  B.  Statistically  stored  dis¬ 
locations  include  closed  loops  and  dipoles  that  do  not  contribute  to  the  total  Burgers  vector  B“  in  (13).  The  second  internal 
state  variable  measures  the  total  density  of  twin  boundaries,  C  =  y/bijj,  where  rjj  is  the  total  area  of  twin  boundaries  mea¬ 
sured  per  unit  volume  in  configuration  B. 


3.2.  Constitutive  assumptions 


Let  p  and  pg  denote  mass  density  of  the  solid  in  current  and  reference  configurations,  respectively,  related  by  pg  =  pj.  Let 
p  =  pf  =  pgj-^  denote  mass  density  in  configuration  B.  The  forthcoming  thermodynamic  analysis  is  conducted  in  B,  the  nat¬ 
ural  configuration  that  serves  as  an  evolving  reference  configuration  for  the  instantaneous  thermoelastic  response  of  the 
crystal  (Eckart,  1948;  Scheidler  and  Wright,  2001;  Clayton  et  al.,  2004a,  2006;  Clayton,  2005).  The  Helmholtz  free  energy 
per  unit  volume  in  B  is  •?  =  pt//,  with  tj/  the  free  energy  per  unit  mass.  The  free  energy  exhibits  the  dependencies 

¥  =  (16) 

where  6  is  the  absolute  temperature.  Electric  polarization  (e.g.,  of  sapphire)  is  not  addressed  in  (16)  since  the  present  study  is 
limited  to  situations  where  electric  fields  are  absent.  Variables  a,  c,  and  C  are,  respectively,  related  to  densities  of  geometri¬ 
cally  necessary  dislocations,  statistically  stored  dislocations,  and  twin  boundaries,  as  discussed  in  Section  3.1.  Set  {/l}  in¬ 
cludes  each  of  the  twin  volume  fractions,  and 


4 = I  = I 


(17) 
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is  a  finite  elastic  strain  tensor  associated  with  elastic  deformation  tensor  C^,  with  components  referred  to  an  assumed  extrin¬ 
sic  Cartesian  frame  on  B  with  basis  vectors  not  tangent  to  possibly  anholonomic  intermediate  material  lines  (Clayton  et  al., 
2004b,  2005).  In  agreement  with  physical  arguments  of  Sections  2.2  and  2.3,  the  free  energy  does  not  depend  explicitly  on 
plastic  deformation  or  twinning  shears. 


3.3.  Thermodynamics 

The  standard  local  forms  of  the  balance  of  energy  and  dissipation  inequality,  respectively,  each  referred  to  the  reference 
configuration,  are  (Eringen,  1962;  Marsden  and  Hughes,  1983) 


Poe  =  i::E-Vi 


■Q.  +  Por,  L:E-po(^i/r  +  lj0'j-^Vo0-(l)SrO. 


(18) 


In  (18),  e  =  i/^  +  0)7  is  internal  energy  per  unit  mass,  with  tj  the  entropy  per  unit  mass.  Symmetric  second  Piola-Kirchhoff 
stress  Z  is  related  to  first  Piola-Kirchhoff  stress  P  and  symmetric  Cauchy  stress  a  by 


17— ly^noB  rp— 1/^ 


(19) 


Symmetric  tensor  E  =  (1  /2)(F’^F  -  G)  is  the  right  Cauchy-Green  strain,  Vo  is  the  covariant  derivative  on  Bo.  Q.is  the  heat 
flux,  and  scalar  r  denotes  other  heat  sources  of  energy  per  unit  mass,  e.g.,  radiation.  Stress  power  per  unit  intermediate  vol¬ 
ume  can  be  written 


rt-iocre  p 


(20) 


where  Mandel’s  stress  (Mandel,  1974)  is  M  =  C  E,  the  symmetric  elastic  second  Piola-Kirchhoff  stress  is  E  =  J  f  o-F  ,  and 
the  velocity  gradient  pulled  back  to  B  is 

L  =  F®-’LF®  =  F^-’F"^  -f  L'  +  L’’  +  (1  /3)j-'l. 

The  time  rate  of  free  energy  change  per  unit  intermediate  configuration  volume  is 

d 


and  following  from  (21)  and  the  symmetry  of  E  and  E®, 

MfL^^  =  +  Mil’;  +  MtV;  +  (l/3)I/-’M,f . 

Expanding  the  rate  of  W  of  (16)  using  the  chain  rule  (Coleman  and  Curtin,  1967), 


•P  =  ^  :  E®  - 
9E® 


dW  ■  dW  ,  dW  ■  dW  ■  dW^i 


with  summation  implied  over  w  twin  fractions  j,  the  entropy  inequality  in  (18)  can  be  written 


^  dW 


.-^):E^-(N  + 


d‘P\  ■  ^ 


dW\ 


J 


3; 


d¥ 


dW  ,  dW  ,  dW  : 


1, 


(^  +  n-(l'+L  +-^1]  -—  :  a  +  +  -—  (  +  P  +  —  V0  ■  q, 


da 


di 


dc "  dp-’  e 


(21) 

(22) 

(23) 

(24) 

(25) 


where  N  =  pi]  is  the  entropy  per  unit  intermediate  volume,  V0  =  VoOF^'  is  the  intermediate  temperature  gradient, 
q=J  ’FQ  is  the  intermediate  heat  flux,  and 


n  =  M  -  ipi 


(26) 


is  a  (negative)  version  of  Eshelby’s  energy-momentum  tensor  or  Eshelby’s  stress  tensor  (Eshelby,  1975;  Maugin,  1994;  Clay¬ 
ton  et  al.,  2004a,  2006)  mapped  to  configuration  6.  Eollowing  standard  arguments  (Coleman  and  Gurtin,  1967;  Scheidler  and 
Wright,  2001;  Clayton,  2005),  stress-elastic  strain  and  entropy-temperature  relations  are  deduced; 


dW 


N=  - 


dW 


aE“ '  de 

The  rightmost  term  of  (25)  contributes  positively  to  dissipation  upon  prescription  of  the  conduction  law 
q  =  -KVe,  -V0-q  =  V0-KVe  5=  0, 

where  K  is  a  symmetric  and  positive  definite  matrix  of  thermal  conductivity.  Applying  (27)  and  (28), 


(27) 


(28) 


n  :  L'  +  n  :  L 


(29) 
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is  the  reduced  dissipation  inequality.  In  the  absence  of  temperature  gradients,  (29)  requires  that  the  energy  dissipated  by 
twinning,  siip,  and  residual  volume  changes  exceeds  the  rate  of  energy  storage  associated  with  defects,  specifically  geomet¬ 
rically  necessary  and  statistically  stored  dislocations  and  twin  boundaries.  From  (9)-(ll),  energies  dissipated  from  twinning 
and  slip,  respectively,  can  be  written 


VV  II  VV  /  It 

n  :  L'  =  5]  ¥/V,  n  :  L-  =  (1  -  fr)  +  E  E  ^jyj 


i=i 


j=i 


(30) 


where  the  driving  forces  are  resolved  stresses  on  each  twin  or  siip  piane,  acting  in  the  direction  of  shear: 

Specific  heat  at  constant  elastic  strain,  measured  per  unit  volume  in  configuration  B,  is  introduced  as 


_  dE 

c  =  =  -9- 

de^ 


(32) 


where  £  =  pe  is  internal  energy  per  intermediate  volume.  Multiplying  the  first  of  (18)  by and  using  (19)-(24),  (26)-(28), 
(32),  the  energy  balance  can  be  written  in  the  intermediate  configuration  as 


temperature  change 


n:(L'  +  L‘’j  H 

dissipation  from  slip  and  twinning 


—  d^\  I 


dissipation  from  defect  generation 


d'F 

dk  ^dOdi 


di  ^dBdi ) ^ 


ep-.E^  +  v_i^ 

thermoelastic  coupling  heat  conduction 


dc  ^oedc  r 


dp  aeop  y 


strain  energy  of  geometrically  necessary  dislocations  strain  energy  of  statistically  stored  dislocations  surface  energy  of  twin  boundaries  energy  of  lattice  reorientation  from  twinning 


(33) 


where  source  r  of  (18)  is  assumed  absent  in  (33)  and  hereafter.  Thermal  stress  coefficients  in  (33)  are 


and  the  anholonomic  covariant  derivative  in  (33)  is  Va  =  Va  +y'E^fP7 with  Va  =  Note  Va 

compatibility  condition  =  0  hoids  and  VoAidj/dF^)  =  =  0- 


(34) 

Va  only  when 


3.4.  Representative  free  energy  potential 


A  particular  form  of  (16)  is  posited  for  anisotropic  crystais  that  may  undergo  iarge  eiastic  deformations,  temperature 
changes,  twinning,  and  disiocation  accumulation.  The  free  energy  is  decomposed  as 


W=W^E^,B.{f>})+Y{e)  +  C,  9) , 


(35) 


where  accounts  for  the  thermoelastic  response,  Y  accounts  for  the  specific  heat  content,  and  W'^  accounts  for  residual  free 
energy  of  lattice  defects.  The  thermoeiastic  energy  consists  of  three  terms: 


1 


1 


=  2  -  IP%{9  -  9o), 


(36) 


with  the  first  term  in  (36)  accounting  for  materially  linear,  but  geometricaliy  noniinear,  mechanical  effects,  the  second 
accounting  for  materialiy  noniinear  elastic  effects  important  at  high  pressures  (Graham  and  Brooks,  1971;  Thurston, 
1 974),  and  the  third  accounting  for  thermoeiastic  coupling.  Here,  9o  is  a  constant  temperature  at  which  the  iattice  param¬ 
eters  exhibit  their  reference  iengths,  and  remaining  coefficients  in  (36)  are  partiai  derivatives  of  free  energy  at  nuli  eiastic 
strain: 


g2,pE 


9ElisdEli 


60O!ftf5£(/i  _ 


,  r'  =  - 


g2,pE 


d9dE 


(37) 


Coefficients  in  (37)  may  depend  on  temperature;  when  measured  at  a  particuiar  temperature,  these  are  referred  to  as  iso¬ 
thermal  elastic  constants.  Superscripts  4  and  6  denote  fourth-  and  sixth  rank  tensors  of  eiastic  coefficients,  usually  referred 
to  (at  fixed  temperature)  as  second-  and  third-order  isothermal  elastic  constants,  respectively  (Brugger,  1964;  Thurston, 
1974).  Typically,  small  differences  between  isothermal  elastic  coefficients  and  their  isentropic  counterparts  measured  in 
sound  speed  experiments  emerge  from  thermal  expansion  (Thurston,  1974). 

In  anisotropic  soiids,  coefficients  (37)  depend  on  the  orientation  of  the  Bravais  lattice  in  configuration  B.  When  twinning 
takes  place,  orientations  of  the  original  reference  lattice  (parent)  and  each  twin  differ.  Here,  a  straightforward  averaging 
method  is  used  to  define  the  effective  coefficients  for  a  voiume  element  consisting  of  the  parent  and  one  or  more  twins. 
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It  is  assumed  that  elastic  deformation  and  elastic  strain  act  uniformly  over  the  parent  and  twins  comprising  this  vol¬ 
ume  element.  Energy  (36)  is  thus  partitioned  into  contributions  from  the  parent  and  each  twin: 


-h)  +  -fj)  -  pfEl^{B  -  eo)(l  -h) 


j=l  L- 


2 


cf^^E\,  +  -  iE;i%^(e  -  eo) 


P, 


(38) 


where  [-gfer  to  coefficients  for  the  parent  lattice,  and  where  for  each  twin  j, 


(39) 


Stress-strain-temperature  relations  following  from  (27),  (36),  and  (38)  are 


2;a/f  ^  4£«ft^££  ^  ££  _ 


^C^^Et  +  l<^cT^^^-‘^El,E%  -  pf{e  -  00) 


(l-/r)  +  E 

J=1  L 


w  r  1 


E\,El,-i^>\e-eo) 


P, 


(40) 


implying  that  elastic  stress  E  for  a  heterogeneous  (twinned)  crystal  Is  equivalent  to  the  volume  average  of  local  stresses  sup¬ 
ported  by  the  parent  and  each  twin.  Effective  material  coefficients  are  thus 


4£a«5  ^  ''Cg'’^^(l  -fr)  + 


J=1 


J=1 


(41) 


'  =  C(1  -/,)  +  J2l}fp,  =  C(1  -/r)  +  J2^(fP- 

j=i  j=i 


From  (41),  necessity  of  inclusion  of  twin  fractions  in  free  energy  functions  (16)  and  (35)  is  evident,  since  effective  ther¬ 
moelastic  coefficients  depend  on  evolving  twin  fractions.  For  anisotropic  crystals  undergoing  twinning,  the  effective  matrix 
of  thermal  conductivity  coefficients  in  (28)  and  (33)  can  be  approximated  in  the  same  manner  as  as  indicated  in  the  last 
of  (41 ),  where  KJ'‘  =  K^^Qj.'^Cps  is  the  conductivity  tensor  of  reoriented  twin  fraction  j  and  is  the  conductivity  tensor  of  the 
parent.  Possible  influences  of  dislocation  densities  (a  and  {)  and  twin  boundaries  (f)  on  effective  elastic  moduli  that  may 
arise  in  some  crystals  at  large  deformations  and  large  defect  densities  (Smith,  1953;  Chung  and  Clayton,  2007)  are  precluded 
by  (35)  since  couplings  between  elastic  strain  and  defects  (oi,  are  not  included.  Such  effects  could  be  Incorporated  by 
generalization  of  (35)  if  deemed  relevant  (Clayton  et  al.,  2004a,b).  The  associated  rate  of  thermoelastic  free  energy  change 
from  rates  of  twin  fractions  is 


dW^  ■■  - 


4  _1fE  /4; 
~  2 


I  ^  f£ 


-  ^cT^^*)e\,e%  -  (/If  -  pf)E%{e  -  Bo), 


(42) 


with  summation  applied  over  j.  Consider  a  situation  in  which  strains  =  0Ccip(B  -  Bo)  arise  from  temperature  change.  The 
following  relationship  emerges  between  the  thermal  stress,  thermal  expansion  constants  («„;;),  and  elasticity  coefficients: 


-  0o), 


(43) 


where  the  second  term  on  the  right  is  negligible  for  most  crystals  with  small  coefficients  <  1-  Symmetry  conditions 
4c”W'5  =  4(Q(a/)){xi5)  ^  6£(aH(x5)(E«_  ^a/)  ^  ^(a/!)  follow  automatically  from  (37).  Voigt’s  notation  (Brugger,  1964;  Thur¬ 

ston,  1974)  exploits  symmetry  of  the  elastic  stress  E,  elastic  strain  E^,  and  these  elastic  coefficients,  where  pairs  of  indices 
11  ^  1,22  ^  2.33  ^  3.23  ^  4,13  ^  5,  and  12  ^  6,  and  1^1'  ^  Ef  1  +  0a/(),  where  a,j5  =  1,2,3  and 

A  =  1 , 2, ...  6.  Relations  between  stress,  strain,  and  temperature  change  in  (40)  then  can  be  expressed  compactly  as 


Ea  =  CabEI  +  2  C/acEfiEc  “  ~  ^o), 


(44) 


where  Cab  and  Cabc  are  second-  and  third-order  elastic  coefficients  written  in  the  standard  notational  scheme  of  Brugger 
(1964)  and  summation  proceeds  over  duplicate  covariant  indices.  Thermal  energy  in  (35)  is  prescribed  as  (Clayton,  2005, 
2009) 


Y  =  -c0ln(0/0o)- 


(45) 
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Finally,  the  residual  energy  of  lattice  defects  in  (35)  is  specified  as 

:  a)*"  +  ,  (46) 

where  /,(  is  an  elastic  shear  modulus  that  may  depend  on  temperature  and  K\,K2,  ■  ■  -  Kg  are  dimensionless  constants  that 
scale  energies  associated  with  internal  state  variables.  The  proper  choice  of  /r  for  anisotropic  materials  is  discussed  in  Appen¬ 
dix  C.  Also,  /  is  scalar  with  dimensions  of  length,  required  from  dimensional  considerations  in  gradient  theories  (Fleck  et  al., 
1994;  Bammann,  2001;  Reguelro  et  al.,  2002;  Abu  Al-Rub  and  Voyiadjls,  2004;  Clayton  et  al.,  2004a,b),  and  N  is  a  constant. 
Recalling  from  Section  3.1  that  ^  the  first  term  on  the  right  of  (46)  provides  a  linear  dependence  of  residual  energy 

on  the  line  density  of  statistically  stored  dislocations,  following  Bammann  (2001 ),  Regueiro  et  al.  (2002),Clayton  et  al.  (2004a, 
2006),  and  references  therein.  Simple  arguments  lead  to  2Ep  =  ic^pib^,  where  Ep  is  the  total  energy  per  unit  line  length  of 
statistically  stored  dislocations.  Including  self-  and  Interaction  energies,  core  energy,  and  stacking  fault  energy  if  the  dislo¬ 
cations  are  partial.  Recalling  that  C  =  the  second  term  provides  for  a  linear  dependence  of  residual  energy  on  the  area 

per  unit  volume  of  twin  boundaries  Similarly,  2Et  =  K2iJ.b,  where  Ej  Is  the  twin  boundaiy  energy  per  unit  area.  In  many 
crystals,  2Et  ^  Wsf,  where  Wsp  is  the  intrinsic  or  extrinsic  stacking  fault  energy  (HIrth  and  Lothe,  1982;  Bernstein  and  Tad- 
mor,  2004).  The  third  term  accounts  for  the  total  line  and  interaction  energy  of  geometrically  necessary  dislocations,  for 
which  a  number  of  forms  for  /,  N,  and  Ks  have  been  suggested  for  different  materials  and  applications  (Fleck  et  al.,  1994; 
Regueiro  et  al.,  2002;  Abu  Al-Rub  and  Voyladjis,  2004;  Clayton  et  al.,  2004a;  Chung  and  Clayton,  2007).  The  fifth  term  ac¬ 
counts  for  Interaction  energies  between  statistically  stored  and  geometrically  necessary  dislocations.  The  fourth  and  sixth 
terms  in  (46)  reflect  interaction  energies  between  twin  boundaries  and  dislocations,  for  example  energies  of  dislocation  lines 
may  be  amplified  at  the  stress  concentration  caused  by  a  pile-up  at  a  twin  boundary  (Christian  and  Mahajan,  1995). 

Methods  for  determining  content  of  geometrically  necessary  versus  statistically  stored  dislocations  have  been  forwarded 
In  recent  years  (Arsenlis  and  Parks,  1999;  El-Dasher  et  al.,  2003;  Hughes  et  al.,  2003),  which  could  facilitate  unique  selection 
of  parameters  in  (46)  if  defect  energies  are  known.  For  example,  energies  can  be  determined  through  measurements  of  cold 
work  (Clarebrough  et  al.,  1957;  Taherl  et  al.,  2006)  or  atomic  simulations  of  defect  energy  (Chung  and  Clayton,  2007).  How¬ 
ever,  in  many  cases,  only  the  total  scalar  line  density  pj  of  all  dislocations,  both  geometrically  necessary  and  statistically 
stored,  may  be  known  from  historical  data.  In  such  cases,  it  becomes  prudent  to  employ  a  reduced  form  of  (46)  written 
in  terms  of  pp. 

(47) 

with  J  =  b^Pp^.  Comparison  of  (46)  and  (47)  implies 

Pt^Ps  +  Pg^  Pc  =  ^  i  =  N=l/2,  =  K3,  K4  =  Ke,  Ks  =  0.  (48) 

Relations  (48)  provide  for  an  equal  contribution,  per  line  length,  of  geometrically  necessary  and  statistically  stored  dis¬ 
locations  to  the  free  energy.  For  uniform  straight  dislocations  of  density  p'  with  identical  Burgers  vectors,  (15)  and  (48)  con¬ 
sistently  give  P(;  =  f)“'(p'b“d'p'ba5^,)'^^  =  p'b“\b“ba)'^^.  Using  (35)  and  (47),  the  rate  of  temperature  increase  in  (33)  can  be 
written’ 


c'e  =  +  e 


[-M] 


+  V  •  KV0, 


where 


(49) 


8'  =  1  - 


+fC4^)p 


Dt  +  - 


.3/2  Pt 


+ 


A 


P}w^ 


(50) 


is  the  Taylor-Quinney  parameter  (Taylor  and  Quinney,  1934;  Clayton,  2005),  such  that  ^  -  p  is  the  ratio  of  inelastic  stress 
power  =  n  :  [L'  +  ?  +  (1  /3)J'J“’  1]  converted  to  stored  energy.  Definition  (50)  is  framed  in  terms  of  the  total  dislocation 
density.  A  more  specific  form  of  (50)  accounting  for  distinct  contributions  from  geometrically  necessaiy  and  statistically 
stored  dislocations  can  be  determined  in  a  straightforward  manner  by  substituting  (35)  and  (46)  into  (33). 


3.5.  Kinetics 


Rate  dependent,  i.e.,  viscoplastic,  kinetic  laws  are  suggested  to  describe  shearing  from  slip  (Hutchinson,  1976;  Teodosiu 
and  Sidoroff,  1976;  Asaro,  1983)  and  twinning  (Wu  et  al.,  2007); 


7' 


yj  =  ys 


P 


yr 

yJ 


0 


(51) 


^  In  a  previous  paper  (Clayton.  2009),  an  incorrect  factor  of  3  preceded  j  in  the  denominator  of  (49).  That  misprint  is  corrected  here  and  is  consistent  with 
(33),  which  has  also  been  corrected  in  a  similar  manner. 
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In  (51 ),  ys  and  yT-  are  material  parameters  with  dimensions  of  1  /t,  m  and  p  are  dimensionless  parameters,  and  g’  and  g>  are 
evolving  resistances— positive  scalars  with  dimensions  of  stress— to  deformation  in  the  parent  grain  by  slip  on  system  i  and 
twinning  on  system  j,  respectively.  In  the  second  of  (51 ),  g]  denotes  resistance  on  slip  system  i  within  reoriented  twin  j.  In  the 
third  of  (51),  2(tJ)  =  V  +  ItiI.  From  (30)  and  (51),  these  rates  are  always  thermodynamically  dissipative  since  for  each  slip 
system,  z'y'  ^  0  and  Tjyj  >  0,  and  since  for  each  twin,  the  product  pziy  >  0.  Because  p  =  0  for  fi  C  0,  the  unidirectional 
nature  of  twinning  is  respected.  In  the  limit  m  ^  oo  or  p  oo,  rate  independent  behavior  is  attained,  respectively,  for  slip 
or  twinning.  In  a  general  sense,  slip  and  twin  resistances  evolve  as 

g^  ^gi(E\e,i,i,c,{p},{g}), 

where  the  rates  depend  not  only  on  the  set  of  state  variables  that  explicitly  enter  the  free  energy  (16),  but  also  on  the  set  of 
hidden  variables  {g}  =  Ig'.ghgjj  with  i  =  1 , . . .  n  and  j  =  1 , . . .  w,  consisting  of  all  possible  slip  and  twinning  resistances.  Slip 
and  twin  resistances  {g}  are  also  referred  to  as  flow  stresses.  Residual  volume  change  J  =J{9,  a,  i(,  {/J})  depends  on  defect 
content  and  temperature,  the  latter  dependence  a  result  of  temperature  influences  on  the  elastic  coefficients  and  their  pres¬ 
sure  derivatives.  Specific  formulae  quantlfyingj  are  deferred  to  (65)-(67)  in  Section  4.4. 

Evolution  equations  for  scalar  Internal  state  variables  c  and  C,  reflecting  respective  densities  of  statistically  stored  dislo¬ 
cations  and  twin  boundaries,  complete  the  model.  An  explicit  evolution  equation  for  the  geometrically  necessaiy  dislocation 
tensor  is  not  needed  since  a  follows  directly  from  kinematic  definition  (15).  Analogous  to  hardening  in  (51 ),  defect  densities 
may  depend  upon  history  of  slip  and  twin  activity,  as  well  as  crystal  structure  and  material  composition.  Generic  evolution 
equations  are 

t  =  C  =  (53) 

From  (52)  and  (53),  free  energy  (16)  thus  depends  on  hidden  variables  {g}  Implicitly,  via  their  influence  on  evolution  of 
the  state  variables. 

Possible  Impedance  or  facilitation  of  slip  or  twinning  via  slip-slip  interactions,  slip-twin  interactions,  and  twin-twin 
interactions  depends  in  a  complex  manner  on  a  number  of  factors.  Including  geometrical  relationships  between  interacting 
systems,  temperature,  crystal  structure,  and  defect  content  (Christian  and  Mahajan,  1995;  Lu  et  al.,  1998;  Castaing  et  al., 
2002;  Wu  et  al.,  2007).  Experimental  data  enabling  unique  quantification  of  these  effects  is  often  scarce,  and  mechanisms 
responsible  for  hardening  are  not  fully  understood  In  many  materials  (Lagerlof  et  al.,  1994).  Because  of  the  large  number 
of  parameters  required  for  a  complete  description  of  interactions  between  individual  deformation  mechanisms,  many  exper¬ 
iments  may  be  required,  with  delineation  of  effects  of  a  particular  mechanism  difficult.  Initial  values  of  hardness  In  (52)  may 
differ  among  different  slip  and  twin  families  in  a  crystal,  and  may  account  for  periodic  lattice  resistance  in  an  initially  perfect 
crystal  (Pelerls,  1940;  Nabarro,  1947),  friction  stress  (Beltz  et  al.,  1996),  and  other  initial  barriers,  for  example  those  resulting 
from  interstitials  in  crystals  with  Impurities  (Kocks  et  al.,  1975).  Short  range  barriers  to  dislocation  motion  are  often  strongly 
temperature-dependent  (Kuhlmann-Wildsorf,  1960;  Kocks  et  al.,  1975),  since  increases  in  temperature  correlate  with  in¬ 
creased  probability  of  dislocations  overcoming  such  barriers  via  thermal  activation.  Long  range  barriers  typically  arise  from 
Interactions  of  local  stress  fields  between  defects,  and  generally  increase  with  defect  densities  that  tend  to  accumulate  with 
strain  until  saturation.  Twinning  resistance  is  often  affected  by  temperature,  and  at  low  temperatures  may  increase  less  stee¬ 
ply  than  slip  resistance  as  the  temperature  is  decreased  (Christian  and  Mahajan,  1995).  This  phenomenon  explains  the  ten¬ 
dency  for  some  solids  to  favor  twinning  over  slip  at  low  temperatures.  Resistance  to  twinning  may  even  decrease  with 
decreasing  temperature  (Christian  and  Mahajan,  1995).  Recent  models  addressing  slip,  twinning,  and  slip-twin  interactions 
via  detailed  hardening  laws  were  advanced  by  Beyerlein  and  Tome  (2008),  Proust  et  al.  (2009),  Shiekhelsouk  et  al.  (2009). 
Proust  et  al.  (2009)  also  modeled  mechanical  detwinnlng. 

The  present  model  assumes  that  thermoelastic  deformation  (and  hence  elastic  strain)  is  uniform  over  the  parent  and  all 
twins  comprising  a  volume  element  at  a  given  “material  point”,  as  Is  clear  from  developments  in  Section  3.4.  Such  an  ap¬ 
proach,  whereby  the  total  stress  at  a  material  point  is  a  volume  average  of  stresses  in  parent  and  individual  twins,  is  similar 
to  earlier  models  of  Kalidlndi  (1998,  2001).  Advantages  and  disadvantages  of  this  type  of  approach  are  discussed  by  Proust 
et  al.  (2007).  Specifically,  the  present  approach  allows  for  straightforward  application  of  hardening  laws  such  as  (52).  How¬ 
ever,  secondary  twinning  and  detwinnlng  are  not  addressed  in  the  present  model.  Furthermore,  stress  equilibrium  among 
twinned  and  parent  regions  with  different  anisotropic  elastic  constants  is  not  strictly  maintained  according  to  the  present 
model.  Alternative  approaches  accounting  for  various  equilibrium  and  compatibility  constraints  have  been  formulated,  using 
micromechanics  arguments,  for  crystals  with  lamellar  structures  (Lebensohn,  1999;  Proust  et  al.,  2007).  Self-consistent  mod¬ 
els  have  also  been  applied  to  anisotropic  polycrystals  with  hexagonal  (Lebensohn  and  Tome,  1993)  and  trigonal  (Tome  et  al., 
1991b)  symmetry.  It  Is  noted  that  the  present  theory  Is  intended  to  be  applied  at  material  points  associated  with  a  small 
length  scale  within  a  single  ciystal,  and  no  assumptions  are  made  regarding  interactions  among  grains  in  a  polycrystal. 
The  present  model  also  Includes  a  number  of  features  often  absent  In  other  treatments;  thermodynamics  of  heat  conduction 
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or  adiabatic  temperature  rise,  thermal  expansion,  higher-order  elastic  constants,  internal  state  variables  associated  with 
energies  of  lattice  defects  such  as  dislocations  and  stacking  faults,  and  possible  volume  changes  induced  by  defects. 


4.  Application:  corundum  single  crystals  with  impurities 

The  framework  of  Section  3  is  specialized  in  Section  4  to  describe  elasticity  and  inelasticity  in  pure  and  doped  corundum. 
Crystal  structure  and  properties,  yield  mechanisms,  strain  hardening,  defect  accumulation,  and  residual  volume  changes  are 
considered.  Initial  yield  mechanisms  are  addressed  for  pure  alumina  single  crystals  over  a  range  of  temperatures,  and  for 
doped  corundum  at  high  temperatures.  The  quantitative  treatment  of  strain  hardening  and  defect  accumulation  focuses  first 
on  pure  and  doped  corundum  deformed  by  basal  glide  at  high  temperatures.  New  predictions  are  then  provided  for  rhom- 
bohedral  twinning  and  for  basal  double  slip  followed  by  rhombohedral  twinning  in  pure  corundum.  New  comparisons  of 
model  features  to  results  of  atomic  calculations  performed  elsewhere  are  reported. 

Some  discussion  in  Sections  4.1  and  4.2  is  abbreviated  from  a  previous  paper  (Clayton,  2009),  wherein  a  more  complete 
description  of  properties  of  pure  alumina  single  crystals  was  given,  and  other  calculations  of  the  response  of  pure  single  crys¬ 
tals  of  various  orientations  to  impact  loading  were  reported.  Discussion  and  corresponding  tables  and  figures  in  Sections.  4.3 
and  4.4  contain  refinements  to  the  previous  model  of  pure  corundum  (Clayton,  2009),  new  results  for  basal  slip  and  defect- 
induced  volume  changes  in  doped  corundum,  and  new  treatment  of  twin  boundaries.  Comparisons  with  atomic  model  re¬ 
sults  in  Section  4.5  are  also  new. 

4.1.  Crystal  structure  and  properties 

The  atomic  structure  of  pure  a-corundum  is  depicted  in  Fig.  3,  following  Kronberg  (1957).  The  primitive  unit  cell  is  rhom¬ 
bohedral  (two  formula  units  per  cell,  lattice  parameter  Qq  =  0.512  nm,  bond  angle  55.3°).  The  hexagonal  unit  cell,  while  con¬ 
sisting  of  more  atoms  (six  formula  units  per  cell),  is  more  convenient  for  describing  mechanical  behavior.  Two  types  of 
hexagonal  cells  are  encountered  in  the  literature  (Kronberg,  1957).  One  is  the  morphological  cell  shown  in  Fig.  3,  with  lattice 
parameters  Ai  =  0.475  nm  and  Cm  =  0.649  nm.  The  second  is  the  structural  cell,  consisting  of  a  second  stack  of  each  of  three 
layers  of  A1  cations  and  O  anions,  differing  from  the  morphological  unit  cell  by  a  rotation  of  180°  about  the  c-axis  [0001], 
Parameters  of  the  structural  hexagonal  cell  are  =  0.475  nm  and  Cs  =  2Cm  =  1.297  nm,  shown  in  Fig.  4.  Elastic  constants 
are  listed  in  Table  1.  In  Table  1  and  henceforth,  all  anisotropic  properties  are  referred  to  the  structural  unit  cell.  Relevant  bulk 
physical  properties  with  supporting  references  are  listed  in  Table  2. 

Prominent  slip  and  twin  systems  are  listed  in  Tables  3  and  4,  respectively,  and  are  illustrated  in  Fig.  4.  While  basal  and 
prism  slip  can  occur  in  either  direction  for  a  given  system,  pyramidal  slip  is  thought  unidirectional  (Heuer  et  al.,  1998)  and 
only  occurs  when  the  resolved  shear  stress  acts  in  a  positive  sense  with  respect  to  the  c-axis,  e.g.,  tensile  loading  along 
[0001].  Twinning  is  also  unidirectional;  for  example,  rhombohedral  twinning  occurs  only  when  the  resolved  shear  stress 
acts  in  a  negative  sense  with  respect  to  the  c-axis,  such  as  occurs  in  compressive  loading  along  [0001],  Thus,  pure  tensile 
deformation  along  [0001]  can  only  be  accommodated  by  pyramidal  slip,  elasticity,  or  fracture,  while  pure  compression 
along  [0001]  can  only  be  accommodated  by  rhombohedral  twinning,  elasticity,  or  fracture.  Slip  directions  for  unidirectional 


Fig.  3.  Morphological  unit  cell  of  sapphire  (Kronberg,  1957). 
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Fig.  4.  Prominent  slip  and  twin  systems  in  corundum  (structural  unit  cell  notation). 


Table  1 

Elastic  coefficients  at  room  temperature. 


Parameter 

Value  [GPa]® 

Parameter 

Value  [GPa]‘‘ 

Parameter 

Value  [GPa]“ 

C„ 

498 

C„, 

-3780 

Cl34 

131 

^12 

163 

C„2 

-1090 

Cl  44 

-302 

Cl  3 

117 

C„3 

-963 

C,55 

-1160 

Cl4 

23 

Cll4 

-55 

C222 

-4520 

C33 

502 

C,23 

-289 

C333 

-3340 

C44 

147 

Cl24 

39 

C344 

-1090 

C,33 

-922 

C444 

19 

®  Winey  et  al.  (2001). 


mechanisms  are  tabulated  here  such  that  the  resolved  shear  stress  for  that  slip  or  twin  system  must  be  positive  to  enact 
shear  on  that  system.  Bourne  et  al.  (2007)  observed  prism  and  basal  dislocations,  twins,  cleavage  fracture,  and  grain  bound¬ 
ary  fracture  in  specimens  recovered  from  impact  experiments  on  polycrystalline  alumina,  with  activity  or  inactivity  of  cer¬ 
tain  mechanisms  dependent  on  the  impact  stress.  Atomistic  simulations  of  hypervelocity  impact  of  sapphire  have  also 
predicted  basal  and  pyramidal  slip  and  basal  and  rhombohedral  twinning  (Zhang  et  al.,  2007,  2008).  Quantification  of  resis¬ 
tances  for  slip,  twinning,  and  fracture  in  corundum  is  currently  an  area  of  active  research  (Tymiak  and  Gerberich,  2007;  Clay¬ 
ton,  2009;  Tymiak  et  al.,  2009).  Because  of  alumina’s  brittleness  at  low  temperatures,  experimental  measurements  of  yield 
mechanisms  at  low  temperatures  must  occur  at  high  pressures  that  suppress  tensile  fracture,  e.g.,  indentation  (Tymiak  and 


1372 


J.D.  Clayton/ International  Journal  of  Plasticity  26  (2010)  1357-1386 


Table  2 

Bulk  physical  properties  (atmospheric  pressure). 


Parameter 

Value 

Remarks 

Po 

3980  kg/m^ 

Room  temperature  mass  density^ 

c 

780 +  0.3(0- 300)  J/kgK 

Specific  heat  per  unit  mass'’ 

an 

5.0  X  lO-'^/K 

Thermal  expansion  coefficient^ 

^33 

5.7  X  10-®/K 

Thermal  expansion  coefficient^ 

Om 

2325  I< 

Melting  temperature'’ 

K 

30  W/mK 

Thermal  conductivity‘s  (300  K) 

P 

156-7O(0/flM)GPa 

Rhombohedral  shear  modulus  (C44)‘* 

^  Burghartz  and  Schulz  (1994). 

Castanet  (1984). 

'  Cahill  et  al.  (1998). 

Zouboulis  and  Grimsditch  (1991). 


Table  3 

Slip  systems  and  slip  dislocations. 


Type 

Direction 

Plane 

Burgers  vector^ 

h  [nm] 

Remarks 

Basal 

[112  0] 

(0001) 

1/3[1120] 

0.475 

Bidirectional'’ 

[1210] 

(0001) 

[2110] 

(0001) 

Prism 

[1100] 

(112  0) 

1/3[1100] 

0.274 

Bidirectional'’ 

[lOlO] 

(1210) 

[0110] 

(2110) 

Pyramidal 

[1101] 

(1011) 

1/3[T101] 

0.512 

Unidirectional 

[0111] 

(1011) 

[1011] 

(0111) 

[1101] 

(0111) 

[0111] 

(1101) 

[1011] 

(1101) 

®  Kronberg  (1957),  Tressler  and  Barber  (1974),  and  Lagerlof  et  al.  (1994). 

Kronberg  (1957),  Snow  and  Heuer  (1973),  and  Heuer  et  al.  (1998). 

'  Tressler  and  Barber  (1974). 


Table  4 

Twin  systems,  partial  dislocations,  and  stacking  fault  energies. 


Type 

Direction 

Plane 

Shear 

Burgers  vector^ 

b  [nm] 

WsF  [)/m2[ 

Basal  (type  ll'’) 

[1010] 

(0001) 

0.635 

1/3[T010] 

0.274 

1. 2-6.7" 

[liOO] 

(0001) 

[0110] 

(0001) 

Rhombohedral  (type  T) 

[foil] 

(1012) 

0.202 

1/21.9[1101] 

0.071 

0.25^ 

[0111] 

(0112) 

[1101] 

(1102) 

®  Kronberg  (1957),  Geipel  et  al.  (1994),  and  Heuer  et  al.  (1998). 

Veit  (1921)  and  Kronberg  (1957). 

'  Heuer  (1966). 

'*  Kenway  (1993)  and  Marinopoulos  and  Elsasser  (2001). 

'  Lagerlof  et  al.  (1984). 


Gerberich,  2007;  Tymiak  et  al.,  2009)  or  confined  compression  (Graham  and  Brooks,  1971;  Castaing  et  al.,  1981;  Scott  and 
Orr,  1983;  Lankford  et  al.,  1998). 

4.2.  Yield  in  corundum 

Flow  stresses  {g}  entering  (51)-(53)  are  decomposed  into  sums  of  contributions  of  various  mechanisms,  e.g.,  following 
Rocks  et  al.  (1975)  and  Clayton  (2005,  2009); 

g'^g‘s+g‘L,  gi^gis+giu  g’-g{+g{- 


(54) 
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Here,  reflects  the  initial  yield  stress  for  slip  in  the  parent  crystal  on  system  i,  and  g[  reflects  long-range  barriers 

associated  with  defects  that  accumulate  during  the  deformation  history.  Analogous  definitions  apply  for  resistance  to  slip 
within  the  twins  and  resistance  to  twin  propagation  in  the  second  and  third  of  (54),  respectively.  Both  terms  in  each  sum 
in  (54)  depend  on  temperature:  the  first,  i.e.,  initial  yield,  is  reviewed  here  for  commercially  pure  corundum  over  a  range 
of  temperatures,  following  Clayton  (2009).  Such  a  review  is  necessary  to  provide  sufficient  background  for  the  discussion 
of  new  results  given  in  Sections.  4.3  and  4.4.  The  second  term  on  the  right  side  of  the  first  equality  in  (54),  g[,  is  newly  ad¬ 
dressed  in  Section  4.3  for  pure  and  doped  alumina  deformed  plastically  at  high  temperatures. 

The  initial  yield  stress  depends  on  short  range  barriers  such  as  strong  PeieiTs  barriers  in  crystals  with  low  initial  defect 
densities  and  non-metallic  bonds  (Friedel,  1964;  Farber  et  al.,  1993),  and  at  high  rates  also  accounts  for  viscous,  phonon, 
and  electron  drag  (Kocks  et  al.,  1975;  Gilman,  1979).  Short  range  barriers  exhibit  a  strong  dependence  on  temperature  for 
both  slip  and  twinning  in  sapphire  (Lagerlof  et  al.,  1994),  reflected  appropriately  by 

g' =g|,exp(|-/l'^j,  g^=g|,exp(|-;/^j,  (55) 

where  gl  is  an  athermal  yield  stress,  Om  is  the  melting  temperature,  and  is  a  dimensionless  parameter.  Analogous  defini¬ 
tions  apply  for  quantities  associated  with  twin  system)  in  the  second  of  (55).  Appropriate  forms  for  the  temperature  depen¬ 
dence  of  the  second  of  (54),  i.e.,  for  slip  systems  within  twinned  regions,  have  not  been  verified  experimentally  for  sapphire, 
but  it  may  be  reasonable  to  assume  that  short  range  barriers  for  a  slip  system  in  the  parent  are  transferred  to  its  rotated 
counterpart  in  a  twinned  region.  Relation  (55)  does  not  conform  to  the  usual  Arrhenius  form  for  thermally  activated  kinetics 
(Kocks  et  al.,  1975)  but  can  be  justified  physically  in  terms  of  activation  volumes  for  cross-slip  in  alumina  (Lagerlof  et  al., 
1994).  Values  of  gg  and  A  are  compiled  in  Table  5  for  pure  alumina  crystals  for  mechanisms  of  basal,  prism,  and  pyramidal 
slip  and  basal  and  rhombohedral  twinning.  Comparisons  to  experimental  data  are  provided  in  Fig.  5  (Clayton,  2009).  Room 
temperature  experimental  values  follow  from  nonlinear  elasticity  calculations  of  driving  forces  for  inelastic  shearing  mech¬ 
anisms  operative  at  the  experimental  HEL  in  single  crystals  of  various  orientations  (Clayton,  2009)  or  from  values  appropri¬ 
ate  for  indentation  (Tymiak  and  Gerberich,  2007).  High  temperature  data  for  basal  and  prism  glide  are  obtained  from 
Lagerlof  et  al.  (1994).  Data  for  pyramidal  slip  are  obtained  from  Tressler  and  Barber  (1974).  For  basal  twinning,  data  are  ob¬ 
tained  from  Castaing  et  al.  (2004),  and  for  rhombohedral  twinning,  from  Scott  and  Orr  (1983).  For  rhombohedral  twinning, 
the  fits  are  most  appropriate  for  300  C  0  C  1300  K,  for  basal  twinning  most  valid  for  300  C  0  ^  1600  K,  and  for  slip  the 
model  validly  spans  the  range  300  ^  0  ^  2100  K.  Trend  lines  extended  to  the  lower  limit  0  =  0  K  are  extrapolations.  Strain 
rate  dependence  of  twinning  is  not  thought  to  be  significant,  based  on  experimental  observations  (Scott  and  Orr,  1983).  Rate 
sensitivities  of  pyramidal  slip  with  m  w  10  (Tressler  and  Barber,  1974)  and  prism  slip  with  m  »  16  (Castaing  et  al.,  1981) 
have  been  observed  for  sapphire  deformed  at  low  rates  and  high  temperatures.  Pressure  and  rate  dependence  are  not  con¬ 
sidered  in  Fig.  5  but  could  be  included  in  a  more  refined  model  when  supporting  experimental  data  spanning  a  wide  range  of 
pressures  and  strain  rates  become  available  for  all  slip  and  twin  families.  Predictions  at  a  temperature  of  300  K  follow  the 
sequence 

pyramidal  slip  basal  twin  basal  slip  prism  slip  _.rhombohedral  twin  /ret 

Ss  ^  ss  ~  Ss  ^  Ss  ^  ss  1 

in  general  agreement  with  trends  for  slip  resistances  in  indentation  (Tymiak  and  Gerberich,  2007). 

4.3.  High  temperature  hardening  and  defect  accumulation  in  pure  and  doped  alumina 

Conceivably,  long-range  barriers  in  each  of  (54)  could  depend  in  a  complex  manner  upon  activity  of  the  slip  or  twin  sys¬ 
tem  under  consideration  (i.e.,  self-hardening)  as  well  as  the  activity  of  all  other  slip  systems  (i.e.,  latent  hardening)  and  twin 
systems  (Christian  and  Mahajan,  1995;  Wu  et  al.,  2007).  Considered  here  are  phenomena  for  which  supporting  experimental 
data  are  available:  hardening  of  basal  slip  and  hardening  of  rhombohedral  twinning  by  total  dislocation  line  accumulation. 


Table  5 

Slip  and  twinning  parameters  for  commercially  pure  alumina. 


Parameter 

Value 

Remarks 

Parameter 

Value 

Remarks 

So  IGPal 

12.7 

Basal  slip 

g„,  IGPal 

53.6 

Basal  slip 

6.6 

Prism  slip 

CO 

10.7 

Basal  slip 

10.0 

Pyramidal  slip 

fio 

8.0  X  10-^ 

Basal  slip 

11.5 

Basal  twin 

9.3  X  10-^ 

Basal  slip 

3.4 

Rhombohedral  twin 

do 

0.54 

Basal  slip 

A 

8.4 

Basal  slip 

0(1 

0.64 

Basal  slip 

5.6 

Prism  slip 

0.58 

Rhombohedral  twin 

3.7 

Pyramidal  slip 

7.7 

Basal  twin 

10.2 

Rhombohedral  twin 
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Fig.  5.  Shear  stresses  for  initial  yield  by  slip  or  twinning  in  pure  corundum  (Clayton,  2009). 


Addressing  first  the  hardening  of  slip  systems,  i.e.,  the  glide  resistance,  the  rate  of  resistance  to  long-range  barriers  in  (54) 
is  prescribed  as 

S[=h-§-y±\f\,  giU  =  0,  (57) 

V  oM/  /c=i 

where  gj^,  is  a  saturation  stress  (Wu  et  ah,  2007)  and  h'  is  a  hardening  modulus.  These  quantities  depend  on  the  thermody¬ 
namic  state  of  the  crystal  as  follows: 

gM  =  ginexp(|-(o'^j,  h'  (58) 

where  gi^,a)',hg,  and  h\  are  material  constants  that  may  differ  among  different  families  of  slip  systems.  These  constants,  as 
well  as  go  of  (55),  may  also  vary  with  the  concentration  of  impurities  or  dopants  (Klassen-Neklyudova  et  ah,  1970;  Pletka 
et  ah,  1982).  For  example,  substitutional  atoms  such  as  Cr,  Ti,  and  Mg  may  be  present  In  natural  (e.g.,  geologic)  crystals, 
or  may  be  added  to  a  synthetic  ceramic  to  adjust  Its  mechanical,  electrical,  and/or  optical  properties  (Castalng  et  ah,  2002). 

Relations  (57)  and  (58)  provide  a  reasonable  fit  to  the  experimental  basal  slip  data  of  Pletka  et  ah  (1977)  for  pure  sapphire 
over  the  range  1673  C  0  C  1993  K,  as  shown  In  Fig.  6(a),  with  model  parameters  listed  in  Table  5.  In  Fig.  6,  y*”  is  the  cumu¬ 
lative  plastic  shear  from  glide  on  a  single  system.  As  shown  in  Fig.  6(b),  the  model  accurately  describes  the  response  of 
corundum  doped  with  cations  of  Cr^*,  Tp^  Ti"'h  or  Mg^*;  atomic  percentages  of  each  dopant  and  corresponding  hardening 
parameters  are  listed  in  Table  6.  For  pure  alumina  and  all  doped  variants,  X  =  8.4  of  Table  5  is  used.  For  all  model  predictions 
and  experimental  results  (Pletka  et  ah,  1982)  shown  in  Fig.  6(b)  except  those  for  Cr^*  doping,  the  temperature  is  0  =  1793  K. 
For  Cr^*-doped  corundum,  predictions  and  experimental  results  (Pletka  et  ah,  1982)  of  Fig.  6(b)  correspond  to  0  =  1773  K. 

As  is  clear  In  Fig.  6(a),  model  predictions  do  not  provide  an  ideal  fit  to  all  experimental  data.  For  example,  dips  in  exper¬ 
imentally  observed  shear  stresses  at  plastic  shear  strains  on  the  order  of  a  few  percent  are  not  captured  by  the  model  of  (57) 
and  (58).  Pletka  et  ah  (1977)  attributed  such  dips  to  an  increase  in  mobile  dislocation  density  associated  with  dislocation 
sources  and  dipole  break-up.  Furthermore,  the  saturation  stress  at  large  plastic  shear  strain  is  not  predicted  accurately  for 
all  temperatures;  e.g.,  the  model  fit  is  excellent  at  1773  K  for  y’’  >  0.20,  while  the  saturation  stress  is  under-predicted  by 
the  model  by  at  least  20%  at  1673  K.  The  law  used  for  temperature  dependence  of  glide  resistance  in  the  present  model 
has  been  selected  as  a  compromise  between  simplicity  (exponential  and  linear  forms  In  (58),  four  parameters)  and  accuracy. 
More  accurate  temperature-dependent  hardening/softening  laws  could  always  be  attempted  at  the  possible  expense  of  addi¬ 
tional  model  complexity  and  added  fitting  parameters. 

Dislocation  accumulation  takes  place  In  conjunction  with  slip  system  hardening  according  to 

0 

g[  =  -  y^),  a'  —  ,  (59) 

VM 

where  pjg  is  the  initial  dislocation  density  and  «{,  and  oc',  are  dimensionless  constants  that  may  vary  among  families  of  slip 
systems.  Following  Pletka  et  ah  (1977),  (59)  extends  the  theory  of  Taylor  (1934)  of  slip  impedance  from  dislocation  interac¬ 
tions  on  parallel  planes  to  conditions  wherein  the  ratio  between  shear  strength  and  dislocation  density  decreases  with 
increasing  temperature.  Relation  (59)  improves  on  a  previous  relation  (Clayton,  2009)  that  did  not  allow  <  Pro  may 
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Fig.  6.  Shear  stresses  for  slip  on  the  basal  plane:  (a)  pure  alumina  over  range  of  temperatures  and  (b)  pure  and  doped  alumina  at  9  «  1800  K. 


Table  6 

Basal  slip  parameters  for  doped  alumina  (kj  1800 1<). 


Dopant 

Atomic 

go  [GPa] 

gm  [GPa] 

^lo 

a 

Cr3+ 

0.20 

14.0 

100 

1.2  X  10-^ 

9.3  X  10-^ 

0.150 

Ti3t 

0.047 

13.5 

55 

9.0  X  10-^ 

9.3  X  10-^ 

0.067 

Ti4t 

0.047 

35.0 

160 

9.5  X  10-^ 

9.3  X  10-^ 

0.15 

0.002 

12.8 

30 

9.0  X  10-^ 

9.3  X  10-^ 

-0.12 

Pletkaetal.  (1982). 


result  from  dislocation  annihilation.  Parameters  are  listed  in  Table  5  for  pure  alumina,  following  calibration  to  data  obtained 
over  the  temperature  range  1673  C  9  C  1993  I<  from  basal  slip  experiments  at  a  rate  of  f  w  ys  =  1.33  X  lO^Vs  (Pletka 
et  al.,  1977).  Values  of  ao  and  listed  in  Table  5  differ  from  those  listed  previously  (Clayton,  2009)  as  a  result  of  new  relation 
(59).  Predictions  are  compared  with  data  (Pletka  et  al.,  1977)  for  pure  sapphire  in  Fig.  7(a).  Predictions  of  the  model  are  com¬ 
pared  with  data  (Pletka  et  al.,  1982)  for  doped  corundum  deformed  by  basal  glide  at  temperatures  0  «  1800  K  in  Fig.  7(b). 
Parameter  a  is  listed  in  Table  6  for  each  doped  variant  of  alumina,  valid  at  6  a  1800  K.  Additional  experiments  at  various 
temperatures  are  required  to  determine  the  temperature  dependence  of  a  in  the  second  of  (59)  for  each  of  the  doped  vari¬ 
ants.  From  Fig.  7,  dislocation  densities  increase  with  plastic  shear  strain  from  initial  value  pjg  =  10®  m“^(Pletka  et  al.,  1977) 
to  saturation  values  in  the  range  lO'^  m^^  <  Pt  <  tn^^  that  increase  with  decreasing  temperature. 

When  all  systems  in  a  given  family  harden  equally  (Taylor,  1934),  for  example  as  in  (57),  relation  (59)  can  be  inverted  and 
used  with  (57)  to  provide  an  evolution  equation  for  the  total  accumulated  dislocation  density  (Clayton,  2005).  Generalization 
of  (57)-(59)  is  required  to  account  for  potentially  more  complex  interactions  between  active  and  inactive  systems  (i.e.,  self¬ 
hardening  versus  latent  hardening),  and  for  interactions  between  statistically  stored  and  geometrically  necessary  disloca¬ 
tions  (Rezvanian  et  al.,  2007).  Both  h'  of  (58)  and  g[  of  (59)  are  scaled  consistently  by  the  rhombohedral  shear  modulus  p 
that  depends  on  temperature  (Table  2  and  Appendix  C). 

Experiments  also  provide  quantitative  evidence  for  hardening  of  rhombohedral  twin  propagation  by  accumulation  of  for¬ 
est  dislocations  (Castaing  et  al.,  2002).  Such  behavior  is  captured  by 

=  =  (60) 
1=1 

where  slip  on  all  systems  can  contribute  to  the  hardness  of  twin  system),  and  hj  is  an  interaction  matrix  relating  the  hard¬ 
ening  rate  of  slip  system  i  to  that  of  twin  fraction).  Setting  hj  =  /ir  for  all  basal  slip  systems  i  and  rhombohedral  twin  systems 
)  gives  g{  =  hr^g'i  with  the  value  of  Hr  listed  in  Table  5,  fit  to  data  of  Castaing  et  al.  (2002)  for  hardening  of  rhombohedral 
twinning  by  basal  slip. 

Rhombohedral  twinning  is  considered  quantitatively  in  Fig.  8(a).  Loading  conditions  in  Fig.  8(a)  correspond  to  uniaxial 
compressive  stress  applied  normal  to  the  basal  plane  (i.e.,  parallel  to  the  c-axis).  For  such  conditions,  the  only  operative 
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Fig.  7.  Dislocation  density  accumulation  from  slip  on  the  basal  plane:  (a)  pure  alumina  over  range  of  temperatures  and  (b)  pure  and  doped  alumina  at 
ffRilSOOK. 


inelastic  deformation  mechanism  is  rhombohedral  twinning,  which  can  occur  simultaneously  on  any  of  the  three  rhombo- 
hedral  twin  systems  of  Table  4,  each  with  a  Schmid  factor  of  0.45  (Scott  and  Orr,  1983).  Representative  experimental  data 
points  in  Fig.  8(a),  corresponding  to  an  axial  compressive  strain  rate  of  e  =  8.7  x  10“®/s,  display  serrations  associated  with 
load  drops  that  correlate  with  fast  twin  nucleation,  propagation,  and  broadening  (Scott  and  Orr,  1983).  In  Fig.  8(a),  e  is  tbe 
nominal  total  axial  strain  and  t  is  tbe  resolved  shear  stress  acting  on  a  twin  system.  Model  predictions  in  Fig.  8(a)  correspond 
to  T  =  gf  =  g{  computed  via  (55).  Because  no  slip  systems  are  active,  (60)  does  not  result  in  any  hardening  associated  with 
forest  dislocations  accumulated  during  slip,  so  the  present  model  provides  a  constant  twin  resistance  at  a  given  temperature. 
The  model  correctly  predicts  a  lower  twin  resistance  at  the  higher  temperature  observed  in  experiments  at  higher  strains 
(e  >  0.12),  but  does  not  predict  tbe  converse  behavior  at  lower  strains.  It  is  noted  that  present  use  of  (60)  with  parameters 
in  Table  5  is  at  the  upper  limit  of  validity,  with  regards  to  temperature,  of  the  fit  shown  in  Fig.  5  for  rhombohedral  twinning 
because  shear  stresses  required  for  rhombohedral  twinning  do  not  decrease  significantly  as  the  temperature  is  increased 
above  1200  I<  (Scott  and  Orr,  1983;  Castaing  et  al.,  2002).  Tbe  model  (at  1173  K)  bounds  most  experimental  data  in 
Fig.  8(a)  from  above,  but  does  not  capture  undulations  in  stress-strain  behavior  associated  in  experiments  with  sequential 
activation  and  fast  propagation  of  individual  twins  within  a  single  ciystal.  Some  such  undulations  are  artifacts  of  the  limited 
velocity  of  the  testing  machine  (Scott  and  Orr,  1983;  Castaing  et  al.,  2002). 

Basal  slip  followed  by  rhombohedral  twinning  is  considered  quantitatively  in  Fig.  8(b).  Boundary  conditions  applicable  to 
Fig.  8(b)  are  explained  by  Castaing  et  al.  (2002);  a  single  crystal  is  compressed  in  uniaxial  stress  applied  along  the  ab  direc¬ 
tion,  at  an  angle  of  57°  from  [0001]  in  the  (1210)  plane,  in  two  stages.  In  the  first  stage,  the  temperature  is  held  at 
9  =  1567  K,  which  favors  basal  slip  (two  systems  each  with  Schmid  factor  of  0.4),  presuming,  as  noted  above,  that  stress 
needed  for  rhombohedral  twinning  does  not  decrease  significantly  above  1200  1<  (Scott  and  Orr,  1983;  Castaing  et  al., 
2002).  The  specimen  is  then  cooled  to  0  =  1169  K.  In  the  second  stage,  the  ciystal  is  reloaded  in  compression  along  the  same 
direction,  which  at  this  lower  temperature,  favors  rhombohedral  twinning  (two  systems  with  a  Schmid  factor  of  0.15)  be¬ 
cause  basal  slip  resistance  increases  significantly  with  decreasing  temperature.  Model  predictions  for  axial  stress  a  (positive 
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Fig.  8.  Resolved  shear  stresses  for  rhombohedral  twinning  in  pure  alumina  (a)  and  axial  compressive  stress  for  basal  double  slip  followed  by  rhombohedral 
twinning  (see  text)  in  pure  alumina  (b). 
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in  compression)  are  compared  with  experimental  results  in  Fig.  8(b).  The  first  loading  stage  activates  double  basal  slip,  and 
simultaneously  leads  to  hardening  of  rhombohedral  twin  systems  via  (60).  After  elastic  reloading,  the  second  stage  activates 
rhombohedral  twinning  on  two  possible  systems.  The  experiment  (Castaing  et  al.,  2002),  performed  at  a  strain  rate  of 
s  =  2  X  10“^/s,  was  apparently  terminated  after  first  occurrence  of  twinning.  The  model  predicts  constant  twin  resistance 
stress  as  twinning  progresses  in  the  second  stage.  From  Fig.  8(b),  the  model  over-predicts  the  saturation  stress  for  basal  dou¬ 
ble  slip  by  10-20%  (recall  that  parameters  used  in  (57)  and  (58)  and  listed  in  Table  5  are  calibrated  independently  to  single 
slip  data),  but  the  rhombohedral  twin  resistance  after  basal  slip  appears  to  be  accurately  modeled. 

Hardening  of  prism  and  pyramidal  slip  and  basal  twinning  systems  by  dislocation  accumulation  have  not  been  quantified 
but  could  be  addressed  by  straightforward  extension  of  (57)-(60).  Potential  twin-twin  interactions,  hardening  of  slip  sys¬ 
tems  within  twins,  and  hardening  of  slip  by  twinning  also  remain  to  be  quantified  for  pure  and  doped  alumina. 

Stored  elastic  energies  of  dislocation  lines  (energy  per  unit  length)  and  twin  boundaries  (energy  per  unit  area)  are  written 


C  I 

471/C 


K\jih 

2 


K2fib 


(61) 


where  in  the  first  of  (61),  K  accounts  for  the  edge  or  screw  character  of  the  dislocation  line  as  well  as  elastic  anisotropy 
(Eshelby,  1949;  Foreman,  1955),  R  is  the  radial  distance  from  the  dislocation  core,  Rc  is  the  cut-off  radius  for  the  dislocation 
core,  and  £  is  a  correction  that  accounts  for  the  core  energy,  line  curvature,  interaction  energies  from  other  defects  and 
boundaries,  and  stacking  faults  associated  with  partial  slip  dislocations.  As  a  first  approximation  (Hull  and  Bacon,  1984;  Heu- 
er  et  al.,  1998),  Ki  =  2.0//C  is  assumed.  The  value  of  K  depends  upon  the  orientation  of  the  dislocation  line  with  respect  to  the 
anisotropic  crystal  lattice  (Foreman,  1955;  Steeds,  1973);  for  large  numbers  of  dislocation  orientations,  it  becomes  prudent 
to  use  isotropic  approximations  K  =  \  for  screw  dislocations  and  IT  =  1  -  v  for  edge  dislocations.  Recall  from  Section  3.4  and 
Table  4  that  VVjf  Ri  0.25  -  10  J/m^  is  a  typical  stacking  fault  energy  associated  with  twin  systems  in  alumina,  from  which  K2 
in  residual  energy  (46)  and  (61)  can  be  estimated  as  Wsf/i/ib).  Kenway,  1993  suggested  that  the  variation  in  stacking  fault 
energy  with  temperature  in  alumina  is  small,  i.e.,  dWsp/dO  w  -1.1  x  10“'' J/(m^K). 

Consider  now  the  residual  elastic  energy  accumulated  in  the  crystal  during  single  slip  at  constant  temperature.  From  (61 ), 
the  cumulative  value  Jl  of  dissipation  fraction  p'  in  (50)  is  approximated  as 


J  g‘dy‘ 


-  l-ib  Pt/K 


j  g'dy' 


(single  slip). 


(62) 


Predicted  dissipation  fraction  [)  is  shown  in  Fig.  9(a)  for  pure  sapphire  deformed  at  various  temperatures,  and  in  Fig.  9(b)  for 
doped  alumina  single  crystals  deformed  at  approximately  1800  K.  In  each  case,  p  is  computed  via  integration  of  (57)  and  use 
of  (59),  (62),  and  material  parameters  listed  in  Tables  5  and  6.  Since  during  experiments  (Pletka  et  al.,  1977,  1982)  disloca¬ 
tions  of  families  with  various  Burgers  vectors  and  line  orientations  were  generated,  in  the  computations  the  representative 
value  used  for  b  is  the  rhombohedral  lattice  parameter  and  used  for  p  the  rhombohedral  shear  modulus.  Shown  are  results 
for  edge  dislocations  and  screw  dislocations;  in  corundum,  energy  storage  ratio  I  -  P  differs  by  a  factor  of  K  =  0.76  between 
these  two  kinds  of  dislocations.  In  the  experiments,  curved  lines  and  loops  were  also  observed,  whose  energies  are  only 
approximated  by  straight  line  formula  in  the  first  of  (61).  At  all  temperatures,  the  stored  energy  appears  fairly  small  relative 
to  the  plastic  work,  since  p  >  0.75,  in  agreement  with  observed  trends  for  many  engineering  metals  (Taylor  and  Quinney, 
1934;  Clayton,  2005). 


Fig.  9.  Predicted  cumulative  heat  dissipation  fraction  from  slip  on  the  basal  plane:  (a)  pure  alumina  over  range  of  temperatures  and  (b)  pure  and  doped 
alumina  at  0  «  1800  K. 
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It  is  noted  that  ( 61 )  and  (62 )  account  for  dislocation  interactions  only  in  an  approximate  way,  via  choice  of  the  empirical  pro¬ 
portionality  constant  K]  =  2.0/K  (Hull  and  Bacon,  1984;  Heuer  et  al.,  1998).  Specifically,  the  energy  per  unit  volume  of  dislo¬ 
cations  is  assumed  proportional  to  their  total  line  length  per  unit  volume,  l.e.,  energy/volume  =  energy/lengthx 
length/volume  =  Ep  x  pj  =  (pb^  /K)  x  Pp,  with  proportionality  constant  £p  assumed  fixed  at  a  given  temperature.  As  demon¬ 
strated  later  In  Section  4.5  and  Table  8,  the  computed  value  of  £p  =  46.3  nj/m  for  full  edge  dislocations  on  the  basal  plane  at  null 
temperature  from  (61)  is  very  close  to  the  value  from  atomic  simulations  (Bodur  et  al.,  2005),  Ep  =  47.4  nJ/m  at  a  dislocation 
density  p-p  ^  m^^  representative  of  the  application  in  Figs.  6, 7  and  9.  A  more  detailed  approach  accounting  for  differences 

in  dislocation  interaction  energies  at  different  dislocation  densities  would  enable  a  possibly  nonlinear  relationship  between 
dislocation  density  and  defect  energy  per  unit  volume. 

Consider  next  the  energy  of  twin  boundaries.  Geometric  arguments  can  be  applied  to  demonstrate  that  the  area  of  twin 
boundaries  per  unit  volume  of  lamellar  twins  Is  approximately 

W 

(63) 

J=1 


where  £  Is  the  mean  thickness  of  a  twin  comprising  fractionj  with  volume  fraction/h  A  proportionality  relationship  between 
characteristic  twin  thickness  and  twin  volume  fraction  was  also  proposed  by  Proust  et  al.  (2007).  For  propagation  of  a  single 
twin  system,  the  analog  of  (62)  is 


=  1  -  WsFpygp 


(single  twin  system), 


(64) 


where  gi  is  assumed  here  only  a  function  of  temperature,  and  the  second  of  (61)  is  used.  Ranges  of  thickness  spanning 
0.1  pm  C  P  C  500  pm  for  basal  and  rhombohedral  twins  In  sapphire  have  been  reported  (Heuer,  1966;  Scott  and  Orr, 
1983;  inkson,  2000;  Castaing  et  al.,  2004).  From  (64),  dissipation  factor  /i  increases  with  increasing  thickness,  twinning  shear, 
and  twinning  stress,  and  decreases  with  increasing  stacking  fault  energy.  Considering  rhombohedral  twins  (twinning  shear 
y  =  0.202),  a  minimum  thickness  of  0.1  pm,  a  maximum  stacking  fault  energy  of  10  J/m^,  and  the  twin  resistance  provided 
by  the  second  of  (55)  and  Table  5,  (64)  gives  a  lower  bound  of  fi  >  0.999  for  0  C  11  C  2100  K.  Hence,  omission  of  the  pure 
twin  boundary  contribution  to  the  ratio  of  stored  to  dissipated  energy  appears  to  be  a  reasonable  assumption  In  crystal  plas¬ 
ticity  modeling  of  alumina.  However,  interaction  energies  between  twins  and  dislocations  could  still  be  tangible  and  may 
merit  consideration  via  nonzero  K4  and  Ks  in  (46)-(48). 


4.4.  Residual  volume  changes 


Residual  volume  change  J  was  introduced  in  Section  3,  following  (3),  (4)  and  (12).  An  estimation  of  this  quantity  from 
nonlinear  elasticity  theory  is  (Zener,  1942;  Seeger  and  Haasen,  1958;  Teodosiu,  1982;  Clayton  and  Bammann,  2009) 


J-1 


A17 

v;;' 


If^-i 

B\dp 


\dP 


(65) 


where  C  and  B  are  Isotropic  shear  and  bulk  moduli,  p  =  -tr<7/3  is  the  Cauchy  pressure,  and  Wo  and  Ws  are  dilatational  and 
deviatoric  strain  energy  densities  per  unit  volume,  respectively,  arising  from  microscopic  elastic  fields  of  defects  within  a 
volume  element  of  crystal.  Volume  change  in  (65)  is  measured  between  reference  and  intermediate  configurations; 
AV  =  V  -  Vo  with  V  the  volume  In  6.  Partitioning  the  elastic  energies  per  unit  dislocation  line  length  In  the  first  of  (61 )  into 
dilatational  and  deviatoric  parts,  (65)  can  be  written  as  (Seeger  and  Haasen,  1958;  Teodosiu,  1982) 


h{;^-f)BpPT  (screw), 

with  Ep  =  pb^(47iKy'  ln(R/Rc)  w  pb^ /K.  Closed  form  expressions  for  residual  elastic  volume  change  exist  only  for  single 
crystals  of  isotropic  or  cubic  symmetry  (Toupin  and  Rivlin,  1960;  Clayton  and  Bammann,  2009),  and  have  been  applied  in 
analysis  of  volume  changes  from  stored  energy  of  cold  work  in  metal  polycrystals  (Wright,  1982).  An  exact  solution  for 
the  residual  elastic  volume  change,  and  possibly  shape  change,  in  a  highly  anisotropic  (e.g.,  trigonal)  Bravais  lattice  would 
require  integration  over  local  elastic  displacement  gradients  induced  by  defects  within  the  body  (Toupin  and  Rivlin, 
1 960;  Clayton  and  Bammann,  2009),  though  simpler  methods  of  approximation  of  volume  changes  based  on  thermodynamic 
arguments  exist  (Zener,  1942;  Holder  and  Granato,  1969).  Here,  following  previous  studies  (Seeger  and  Haasen,  1958),  the 
isotropic  approximation  (66)  is  investigated  as  an  order-of-magnitude  estimate,  consistent  with  the  isotropic  approximation 
of  dislocation  energies  in  (61).  Formal  derivations  and  limitations  of  (65)  and  (66)  are  discussed  by  Clayton  and  Bammann 
(2009). 

Table  7  shows  values  of  the  normalized  volume  change  AV/(Lb^),  where  L  is  the  total  length  of  dislocation  lines  of  pure 
screw  or  pure  edge  character,  and  b  is  the  magnitude  of  the  Burgers  vector.  For  alumina,  effective  isotropic  constants  (Holm 
et  al.,  1999  and  Appendix  C)  and  their  pressure  derivatives  (Sarkar  et  al.,  1996)  are  obtained  from  the  literature.  Values  of  the 


J-l=^  = 
^  Vo 
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normalized  volume  change  In  AI2O3  from  dislocation  lines  are  positive  in  agreement  with  other  solids  (Seeger  and  Haasen, 
1958;  Wright,  1982),  are  very  to  close  those  for  metallic  crystal  Cu  (face  centered  cubic),  and  are  ~  2/3  of  those  for  ionic 
crystal  NaCl  (cubic  rock  salt).  Data  for  residual  volume  changes  and  dislocation  densities  are  available  for  Cu  polycrystals 
(Clarebrough  et  al.,  1957);  the  range  of  values  for  the  volume  change  per  line  length  in  the  rightmost  column  of  Table  7  cor¬ 
responds  to  applied  compressive  strains  ranging  from  0.3  to  0.7  and  dislocation  line  densities  found  from  energy  measure¬ 
ments.  The  theory  thus  provides  reasonably  accurate  results  for  Cu,  within  a  factor  of  two  of  the  experimental  values. 

Shown  in  Fig.  10(a)  are  volume  changes  computed  from  (66)  and  properties  in  Table  7,  for  distributions  of  pure  edge  and 
pure  screw  dislocations  of  line  densities  predicted  by  the  model  in  Fig.  7(a)  for  basal  slip  at  various  temperatures.  Shown  in 
Fig.  10(b)  are  volume  changes  computed  from  (66)  and  the  properties  in  Table  7,  for  dislocations  in  pure  and  doped  alumina 
of  line  densities  predicted  by  the  model  in  Fig.  7(b)  for  basal  slip  at  0  w  1800  K.  At  fixed  temperature  and  composition,  the 
predicted  ratio  of  residual  expansion  in  alumina  for  pure  screw  dislocations  to  that  for  pure  edge  dislocations  is  w  0.6.  As 
demonstrated  by  Clayton  and  Bammann  (2009),  the  volume  change  for  a  ciystal  containing  a  mixture  of  non-interacting 
edge  and  screw  dislocations  would  fall  between  predictions  for  densities  of  pure  edge  and  pure  screw  dislocations.  Volume 
changes  remain  small  as  saturation  levels  of  dislocation  density  observed  in  high  temperature  basal  glide  are  approached; 
AV/Vo  »  10“^  corresponds  to  the  volume  increase  that  would  result  from  thermal  expansion  in  alumina  from  a  temperature 
rise  of  »  0.7  K. 

According  to  the  theory,  volume  increases  on  the  order  of  1%  would  be  achieved  upon  generation  dislocation  densities 
Pj  Ri  10'®  m^^,  corresponding  to  a  dislocation  spacing  on  the  order  of  10  nm.  Conceivably,  very  large  dislocation  densities 
may  be  required  to  enable  extremely  large  plastic  strain  rates  in  shock  physics  experiments  (Rohatgi  and  Vecchio,  2002) 
if  the  upper  bound  on  dislocation  velocity  is  limited  to  the  velocity  of  transverse  acoustic  waves  (Rocks  et  al.,  1975).  Den¬ 
sities  Pj  a  10'®  m^^  have  been  approached  in  such  cases  (Merala  et  al.,  1988;  Rohatgi  and  Vecchio,  2002).  Generation  of  de¬ 
fect  densities  of  this  magnitude  could  affect  the  measured  pressure-volume  response;  for  example,  the  pressure  required  to 
offset  a  1%  volume  increase  would  be  on  the  order  of  0.016  fs  2.5  GPa.  However,  it  is  noted  that  the  approximations  in  (61) 
and  (66)  may  decrease  in  validity  as  dislocation  spacings  decrease  and  interaction  energies  among  dislocations  become 
stronger.  Atomic  simulations  offer  the  possibility  of  more  accurate  predictions  of  effects  of  large  defect  densities  on  volume 
changes,  effective  (tangent)  elastic  moduli,  and  stored  energies  (Dienes,  1952;  Bodur  et  al.,  2005;  Chung  and  Clayton,  2007). 


Table  7 

Effective  isotropic  elastic  properties  and  volume  changes  associated  with  dislocations. 


Crystal 

G  [GPa] 

B  [GPa] 

V 

dC/dp 

dB/dp 

AV/{Lb^)edge 

AV/{Lb^)  screw 

AV/(Lb^)  experiment‘d 

Al203^ 

157 

250 

0.24 

1.7 

4.2 

1.7 

1.0 

- 

NaCC 

15 

23 

0.24 

2.0 

5.8 

2.5 

1.4 

- 

Cu'’ 

47 

152 

0.36 

1.2 

5.6 

1.8 

1.0 

1.6-2.3 

®  Sarkar  et  al  (1996)  and  Holm  et  al  (1999). 

Seeger  and  Haasen  (1958). 

'  Clarebrough  et  al.  (1957). 
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Fig.  10.  Predicted  volume  changes  resulting  from  dislocation  lines  accumulated  during  slip  on  the  basal  plane:  (a)  pure  alumina  over  range  of  temperatures 
and  (b)  pure  and  doped  alumina  at  6  »  1800  K. 
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For  example,  Bodur  et  al.  (2005)  reported  a  decrease  in  energy  per  unit  length  Ep  of  periodically  arranged  edge  dislocation 
dipoles  on  the  basal  plane  by  about  a  factor  of  3,  from  47  4  nj/m  (Table  8)  to  1 6.5  nj/m,  as  the  dislocation  density  is  increased 
from  pj  fs  lO'^  irf  to  pj  fs  lO’®  m^^.  Accounting  for  this  decrease,  the  volume  increase  computed  via  (66)  at  a  large  disloca¬ 
tion  density  of  pj  rs  lo'®  m^^  would  be  reduced  by  a  factor  of  3,  to  about  0.3%,  corresponding  to  a  pressure  change  on  the 
order  of  0.8  GPa.  A  procedure  for  computing  residual  volume  changes  (as  well  as  possible  shape  changes)  associated  with 
defects  in  a  self-equilibrated  lattice  in  terms  of  anharmonic  contributions  to  inter-atomic  potentials  is  outlined  by  Clayton 
and  Bammann  (2009). 

Finally  consider  volume  changes  attributed  to  twin  boundaries  (i.e.,  appropriate  stacking  fault  energies)  bounded  accord¬ 
ing  to  the  second  of  (61)  and  (65)  as 


J_/^ 

2C\dp 


WspiJj  ^ 


2B  \dp 


1 


WsfPt- 


(67) 


Holder  and  Granato  (1969)  suggested  that  the  lower  bound  in  (67)  would  be  appropriate  for  low-angle  grain  boundaries 
constructed  from  arrays  of  screw  dislocations.  Using  values  from  Table  7  for  elastic  constants  and  their  pressure  derivatives, 
and  using  a  maximum  stacking  fault  energy  of  lOJ/m^,  (67)  becomes  3.4  x  C  AV/\/o  ^  6.4  x  10“'';7^,  with  the 

area  per  unit  volume  of  twin  boundaries  t]j  in  dimensions  of  m“\  From  (63),  withp  =  1  and  a  minimum  twin  thickness 
d  =  0.1  pm  (Inkson,  2000),  bounds  for  the  volume  change  are  6.8  x  10  ®  ^  AV/Vo  C  1.3  x  10“^.  Comparing  with  results 
in  Fig.  10,  this  maximum  volume  change  from  the  stored  energy  of  twin  boundaries  is  of  the  same  order  of  magnitude  as 
that  predicted  from  dislocation  lines  accumulated  in  single  crystals  of  corundum  at  high  temperatures  at  shear  strains  on 
the  order  of  0.25.  Not  considered  here  are  possible  effects  of  doping  on  elastic  coefficients,  which  could  be  important  if  dop¬ 
ing  alters  the  concentration  of  vacancies  or  interstitials  (Dienes,  1952). 


4.5.  Comparison  with  atomic  calculations 

In  Table  8,  parameters  or  predictions  of  the  present  model  are  compared  with  results  of  atomic  calculations  published 
elsewhere  or  with  theoretical  models.  The  initial  resistance  to  slip  of  basal  dislocations  of  the  present  model  in  (55)  is  com¬ 
pared  to  values  obtained  from  molecular  dynamics  (MD)  simulations  of  hypervelocity  impact  (Zhang  et  al.,  2007,  2008)  and 
consideration  of  the  generalized  stacking  fault  energy  gradient  (Nishimura  et  al.,  2009)  using  MD  simulations  as  well  as 
quantum  mechanics  (density  functional  theory  or  DFT).  Also  shown  are  values  of  the  theoretical  Peierls-Nabarro  stress  (Pei- 
erls,  1940;  Clayton,  2009)  for  full  edge  and  screw  dislocations  on  the  basal  plane.  The  value  of  initial  slip  resistance  for  the 
present  model  at  null  temperature,  12.7  GPa,  falls  between  Peierls-Nabarro  stresses  for  edge  and  screw  dislocations,  and  is 
15—35%  lower  than  the  atomic  simulation  results.  It  is  noted  that  the  present  model  does  not  incorporate  effects  of  hydro¬ 
static  pressure  on  the  basal  slip  resistance.  Extremely  high  pressures  generated  during  hypervelocity  impact  (Zhang  et  al., 
2008)  or  indentation  (Tymiak  and  Gerberich,  2007)  could  conceivably  cause  an  increase  in  the  apparent  glide  (as  well  as 
twinning)  resistance,  as  has  been  noted  in  MD  simulations  of  other  materials  (Xu  et  al.,  2004).  For  example,  if  the  glide  resis¬ 
tance  is  assumed  proportional  to  the  shear  modulus  (Peierls,  1940;  Nabarro,  1947;  Rocks  et  al.,  1975),  at  a  pressure  of 
100  GPa,  the  shear  modulus  and  hence  the  glide  resistance  would  approximately  double  in  magnitude  over  its  value  at  ambi¬ 
ent  pressure,  if  the  pressure  derivative  of  the  shear  modulus,  1.7  in  Table  7,  sufficiently  describes  the  variation  of  shear  mod¬ 
ulus  up  to  such  pressures.  In  previous  work  (Clayton,  2009),  the  lower  value  of  4.3  GPa  at  room  temperature  was  found  to 
provide  a  reasonable  explanation  for  possible  occurrence  of  basal  slip  in  plate  impact  experiments  on  single  crystals  (Graham 
and  Brooks,  1971)  at  lower  pressures  on  the  order  of  10-20  GPa. 


Table  8 

Comparison  of  model  features  with  atomic  or  theoretical  calculations. 


Quantity 

Present  model 

Atomic/theoretical  calculation 

Value 

Equation 

Value 

Reference 

Basal  glide  resistance 

12.7  GPa(0  K) 

(55) 

15  GPa 

a 

4.3  GPa  (300  I<) 

(55) 

16.3  GPa(MD) 

18.9  GPa  (DFT) 

b 

9.6  GPa  (edge.  0 1<) 

c 

17.9  GPa  (screw.  0  K) 

c 

Basal  dislocation  energy  Ep 

46.3  nJ/m  (edge,  0  K) 

(61) 

47.4  nJ/m  (edge,  0  K) 

d 

35.2  nJ/m  (screw,  0  K) 

(61) 

Basal  stacking  fault  dilatancy  {AV/V)/t]j 

0.028-0.052  nm 

(67) 

0.028-0.096  nm 

e 

Zhang  et  al.  (2007,  2008). 
Nishimura  et  al.  (2009). 

'  Peierls  (1940)  and  Nabarro  (1947). 
'*  Bodur  et  al.  (2005). 

®  Marinopoulos  and  Elsasser  (2001). 
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The  energy  per  unit  length  of  full  dislocations  on  the  basal  plane  obtained  using  (61 )  with  Kr  =  2/K  and  fi  =  C44  is  com¬ 
pared  in  Table  8  with  results  of  MD  simulations  of  full  basal  plane  edge  dislocations  (Bodur  et  al.,  2005).  The  edge  dislocation 
energy  agrees  closely  with  the  atomic  calculation,  the  latter  valid  for  a  dislocation  line  density  of  Ri  lO’^  m^^  (Bodur  et  al., 
2005)  in  the  null  temperature  limit. 

In  the  final  row  of  Table  8,  the  volume  increase  from  basal  plane  stacking  faults  computed  using  the  bounds  in  (67)  is 
compared  with  results  of  molecular  simulations  (Marinopoulos  and  Elsasser,  2001).  Specifically,  the  volume  change  per  den¬ 
sity  of  stacking  faults  is  computed  in  (67)  using  a  representative  stacking  fault  energy  Wjf  =  4.0J/m^,  multiplied  by  a  factor 
of  two  because  here  rjj  represents  the  stacking  fault  density  rather  than  the  twin  boundary  density.  Atomic  simulation  re¬ 
sults  follow  from  consideration  of  computed  expansions  normal  to  various  kinds  of  basal  plane  stacking  faults  ranging  in 
magnitude  from  5.5%  to  18.7%  of  the  rhombohedral  lattice  parameter  Qq,  computed  for  relaxed  structures  at  null  external 
pressure  and  null  temperature  (Marinopoulos  and  Elsasser,  2001).  Basal  plane  stacking  faults  are  associated  with  partial  dis¬ 
locations  of  type  1/3(1  01 0) {0001},  which  only  involve  the  cation  (i.e.,  Al)  sub-lattice,  since  1/3(1010)  is  a  perfect  trans¬ 
lation  vector  in  the  anion  (i.e.,  O)  sub-lattice.  Agreement  between  the  present  model  and  the  quoted  atomic  calculations  is 
reasonable.  Dilatancy  for  various  relaxed  stacking  fault  structures  on  basal  and  prism  planes  was  also  computed  by  Kenway 
(1993). 


5.  Conclusions 

Anisotropic  mechanisms  of  elasticity,  plastic  slip,  and  twinning  have  been  addressed  in  a  crystal  plasticity  model  with 
consistent  thermodynamics.  Einite  deformation  theory  is  required  for  addressing  large  shears  arising  from  slip  and  twinning, 
lattice  reorientation  arising  from  twinning,  and  nonlinear  elastic  effects.  Accordingly,  each  mechanism  has  been  represented 
explicitly  within  a  multiplicative  decomposition  of  the  deformation  gradient.  A  constitutive  framework  with  internal  state 
variables  provides  thermodynamic  relationships  among  state  variables  and  driving  forces  for  time  rates  of  inelastic  deforma¬ 
tions.  Internal  state  variables  consist  of  geometrically  necessary  and  statistically  stored  dislocation  densities,  twin  boundary 
area  densities,  and  twinned  volume  fractions.  Dissipative  inelastic  deformation  rates  consist  of  rates  of  plastic  shear,  rates  of 
twin  fractions,  and  rates  of  volume  changes  attributed  to  defect  content.  The  present  theoretical  model  extends  that  of  pre¬ 
vious  work  (Clayton,  2009)  by  its  explicit  incorporation  of  geometrically  necessary  dislocations  corresponding  to  gradients  of 
inelastic  stretch  and  rotation  and  statistically  stored  dislocations  associated  with  dislocation  dipoles  and  closed  loops.  Kine¬ 
matic  identities  demonstrate  the  contributions  of  slip  and  twinning  to  the  geometrically  necessary  dislocation  density  tensor 
(e.g.,  full,  partial,  and  interfacial  dislocations).  Both  geometrically  necessary  and  statistically  stored  dislocations  are  ad¬ 
dressed  in  the  thermodynamic  framework. 

A  fully  nonlinear,  anisotropic,  crystal  plasticity  model  for  mechanical  and  thermodynamic  behavior  of  pure  and  doped 
corundum  has  been  developed  and  exercised,  accounting  for  lattice  defects  mentioned  above,  i.e.,  dislocation  densities,  twin 
boundary  densities,  and  twinned  volume  fractions.  Specific  inelastic  deformation  modes  include  basal,  prism,  and  pyramidal 
slip  and  basal  and  rhombohedral  twinning.  Hardening  from  dislocation  accumulation  during  basal  slip  has  been  quantified  in 
pure  and  doped  corundum  single  crystals,  and  stored  energies  associated  with  dislocations  and  twin  boundaries  have  been 
estimated.  The  cumulative  ratio  of  stored  energy  from  dislocations  to  dissipated  energy  in  basal  slip  is  predicted  at  less  than 
25%  for  plastic  shears  of  less  than  0.25  and  temperatures  greater  than  1600  K,  regardless  of  dopant  concentration.  Upon  con¬ 
sideration  of  stacking  fault  energies  and  twin  thicknesses  observed  in  sapphire,  stored  energies  associated  exclusively  with 
twin  boundaries  appear  orders  of  magnitude  smaller  than  those  of  dislocation  lines.  Volumetric  expansions  from  nucleation 
of  dislocation  lines  and  twin  boundaries  have  been  estimated  using  formulae  originating  from  nonlinear  elasticity  theory 
(Clayton  and  Bammann,  2009,  and  references  therein).  Such  expansions  are  predicted  small  in  alumina  for  dislocation  line 
densities  observed  in  low  pressure,  high  temperature  experiments  (Pletka  et  al.,  1977, 1982),  and  are  also  predicted  small  for 
twin  boundaries  of  spacing  (Inkson,  2000)  and  energy  (Lagerlof  et  al.,  1984)  observed  experimentally.  New  predictions  for 
rhombohedral  twinning,  including  hardening  by  prior  basal  slip,  are  provided;  these  predictions  offer  a  reasonable  depiction 
of  experimental  results  (Scott  and  Orr,  1983;  Castaing  et  al.,  2002).  Einally,  model  features  pertinent  to  defects  on  the  basal 
plane  are  compared  with  atomic  calculations.  Expansion  associated  with  stacking  faults  predicted  via  the  nonlinear  elastic 
approach  of  the  present  theory  compares  favorably  with  predictions  of  atomic  simulations  (Marinopoulos  and  Elsasser, 
2001). 


Appendix  A 

Isochoric  characteristics  of  twinning  and  slip  noted  in  (12)  are  derived  in  what  follows.  The  Jacobian  determinant  of  the 
deformation  gradient  in  (2)  is  computed  as  (Truesdell  and  Toupin,  1960;  Eringen,  1962;  Marsden  and  Hughes,  1983) 


^  =  \/il  =  g  (Al ) 


where  by  notational  convention,  determinants  of  metric  tensors  g  and  G  are  absorbed  into  permutation  tensors  Eabc  and 
respectively.  It  follows  that  (Ericksen,  1960;  Thurston,  1974) 
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dF,  6 


_  ^  /o  cd  I  c  c^^BArd  rb  ,  o  ^BG^nb  pC  \  _  ^  „  ^CDArC  pd  _  ^  „  ^CDArC  cd  tE 

—  c  i^cdfl^  ^  C‘D  ^  ^dba^  ‘  D‘  B  ^  ^bca^  ^  c)  —  n  ^cda^  ‘  D  —  n^cda^  ^  C‘D^E 


^  P  qCDAhC  pd  pe  C-IE  ^  „  f,CDET<  pd  pe  p-lA 

^c^cda£  f.c^D'^E^e  =  c^cdeS  ^cF-oF-EFa  ^  JF  a  ■ 


(A2) 


To  verify  the  final  two  equalities  in  (A2),  notice  that  =  3y  =  (1  The  material  time  derivative  of 

(A1 )  is  computed  as  follows  using  the  chain  rule,  noting  that  material  time  derivatives  of  metric  tensors  vanish  for  stationary 
coordinate  systems  (Eringen,  1962): 


i  =  §aK -JF::''F^a -iF^F::''  =7trL, 


where  L  is  the  spatial  velocity  gradient  of  (8).  Applying  the  same  arguments  to  /  of  (4), 

9J'  _  r/p;-l/! 

■la  J  ^  a  ’ 


dF; 

'll  _  Jj!_pl<i  —  [Ipl-lfipl'x  _  jlplapl-lp  _  jljla  _  jlfyfl 
J  ^  r^rla'^  f  ^ -d  a.  -/I  ^  J  a  —  J'-ix— JOL, 
oF.i, 


(A3) 

(A4) 

(A5) 

(A6) 


with  L'  defined  in  (9)  and  (10).  From  (5)  and  (10), 

(W  \  W 

^/A/si,  (g)  mi,  =/  ^py  si,  •  mi,  =  0. 

J-'  /  J=' 

Since  by  definition, /'(t  =  0)  =  1,  (A7)  yields  the  isochoric  conditions 
/  =  l(Vt  55  0). 

Now  consider  the  product  F  =  F'f’’.  The  following  equalities  apply: 

J =j'f  =  Uani^^hH% 


dh 


= 

OL  5 


7  =  =  JFp'^F^  =  JF^aFP''  =7tr(F'F'-’  +  f'FPF'’-’F'-’)  =  //tr(L‘  +  LP), 

dF^A 

where  definitions  in  (9)  have  been  used.  Applying  the  product  rule  to  (A9)  and  using  (A7)  and  (A8)  in  (All)  gives 

7=  7  /+  7  /(trL'+trL'’)=/trLP. 


(A7) 

(A8) 

(A9) 

(AlO) 

(All) 

(A12) 


=0 


=1 


=1 


=0 


Now  taking  the  trace  of  (11)  and  using  (5)  and  (6),  (A12)  reduces  to 

w  /  n 

(’’si,  ®  mi,  + 

1=1 


f  ^/trLP  =/tr 


( 1  -  /t)  E  7®!,  ®  mi,  +  E  7  E  yjsi,,  ®  m'^. 

■■  ’  7=1  V  1=1 


-f 


n  w  /  n 

(1  -/r)  E  7*0  ■  mi,  +  E  [P E  •  >"oj 

i=l  7=1  \  i=l 


(A13) 


/ 


=77i-/r)E7s^+7E 

=0  ^ 


\ 


7E7(Q^si,).(mi,(iA) 

V  ' / 


=  0. 


Since  by  definition, 7(f  =  0)  =  1,  (A13)  yields 
/  =  l(Vf  5=  0). 


(A14) 


j.D.  Clayton/ International  Journal  of  Plasticity  26  (2010)  1357-1386 


1383 


Appendix  B 


An  alternative  to  multiplicative  decomposition  (3)  in  which  the  sequence  of  F'  and  F"*  is  reversed  is  considered  in  what 
follows,  it  is  demonstrated  that  the  sequence  of  F'  and  F**  in  the  multiplicative  decomposition  does  not  significantly  affect 
other  aspects  of  the  theoretical  model  of  Section  3. 

Decomposition  (3)  is  replaced  with 

F  =  F^J’/^FV  =  F'^F  =  F^-F,  F‘•  =  F®J’''^  F  =  F'’f'.  (B1) 


Jacobian  determinants  of  deformation  mappings  in  (4)  are  instead 

J=V^detF  =/]//,  /  =  yd^detF^  /  =  detF‘’,  /  =  \/g^detF'.  (B2) 

Spatial  velocity  gradient  (8)  is  replaced  with 

L  =  Ff^’  =  F^f*^-’  +  f'^F'’f‘’-’f'^-’  +  f'^f'’F'f'-’f'’-’f®-’  +  (1/3)]/-' 1. 

Inelastic  velocity  gradient  (9)  is  instead 

FF-’  =F‘’F'’-Vf'’F'F'-’F‘’-’  +(l/3)i/-’l  ^L'’  +  L'  +  (1/3)JJ-H, 
where  the  following  definitions  apply  instead  of  (10)  and  (11): 

W 

P  =  f‘’F'F'-’f‘’-’  =  ®  mi). 

J=1 

n  w  /  n  \ 

L'’  =  (1  -It)  ®  mj,  +  ^  ^  yjsoj  ®  mlij  . 

i=i  j=i  \  i=i  / 


(B3) 

(B4) 

(B5) 

(B6) 


As  a  consequence  of  (B5)  and  (B6),  and  similarly  to  (12),  twinning  and  slip  are  always  isochoric: 

trr  =  trL'’  =  0^]'=/  =  0^/=/  =  l  (t  Ss  0).  (B7) 

so  that  the  first  of  (B2)  reduces  to  J  =  f]  =  /.  The  total  Burgers  vector  associated  with  geometrically  necessary  defects  in  (13) 
is  replaced  with 


(B8) 


where  the  two-point  dislocation  density  tensor  of  (14)  is  instead 


_  pABCrUic^'y.  ,  pABCpPap^/^ 

slip  gradients  twin  gradients 


(B9) 


The  geometrically  necessary  dislocation  tensor  *  following  the  second  and  third  equalities  of  ( 1 5 )  is  unchanged.  The  treat¬ 
ment  of  constitutive  behavior,  thermodynamics,  and  kinetics  in  Sections  3. 2-3. 5  is  unchanged  except  that  L‘  is  replaced  with 
L'  of  (B5)  and  U  is  replaced  with  P  of  (B6).  Notice  that  the  right  sides  of  (B5)  and  (B6)  are,  by  assumption,  identical  to  right 
sides  of  (10)  and  (1 1 ),  respectively,  so  that  FF^’  is  computed  identically  in  terms  of  slip  and  twin  system  activity  regardless  of 
the  sequence  of  F'  and  F**  in  the  decomposition  of  F. 


Appendix  C 

The  shear  elastic  constant  ji  entering  free  energy  contribution  (46)  can  be  prescribed  in  a  number  of  ways  for  anisotropic 
materials.  When  the  energy  of  lattice  defects  is  dominated  by  a  single  prominent  family  of  dislocations,  /,(  should  be  chosen 
consistently  with  JCi  to  yield  the  correct  energy  per  unit  length  of  dislocation  line.  Various  solutions  for  the  energy  factor  for 
dislocations  in  anisotropic  crystals  are  given  by  Foreman  (1955),  Steeds  (1973).  When  different  kinds  of  dislocations, 
stacking  faults,  and  twin  boundaries  affect  the  energy  as  in  (46),  is  simply  viewed  as  an  effective  (isotropic)  elastic  shear 
modulus,  and  coefficients  Ka,K2,...Ks  then  scale  the  energy  associated  with  each  kind  of  defect  in  an  appropriate  manner. 

A  number  of  averaging  or  self-consistent  methods  are  available  for  computation  of  effective  elastic  constants  (Hill,  1952; 
Sarkar  et  al.,  1996).  Perhaps  the  most  common  is  the  Voigt  average  (Hill,  1952): 

Mv  ~  Yg  [(Cii  +  C22  +  C33)  —  (C12  +  C23  +  C13)  +  3  (C44  +  C55  +  Ces)] ,  (Cl) 

By  =  g  [Cii  +  C22  +  C33  +  2(Ci2  +  C23  +  C13)], 


(C2) 
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where  fty  and  By  are,  respectively,  Voigt-averaged  shear  and  bulk  moduli  that  provide  an  upper  bound  to  the  true  effective 
elastic  constants  of  an  isotropic  polycrystal.  For  sapphire,  noting  that  for  a  trigonal  crystal  Cn  =  C22,  C13  =  C23,  C44  =  C55,  and 
Css  =  (1  /2){C]]  -  C12)  (Thurston,  1974),  (Cl)  and  (C2)  yield  jiy  =  166  GPa  and  By  =  255  GPa.  Various  experimental  and  the¬ 
oretically  computed  values  of  isotropic  shear  and  bulk  moduli  quoted  by  Holm  et  al.  (1999)  for  a-alumina  fall  in  the  range 
156  GPa  C  G  C  166  GPa  and  245  GPa  C  B  C  255  GPa,  with  effective  Poisson’s  ratio  0.23  C  v  ^  0.24. 

For  calculations  in  Sections  4.2  and  4.3  of  the  present  paper  that  focus  primarily  on  basal  slip,  ft  =  C44  has  been  used,  since 
this  is  considered  an  appropriate  shear  elastic  coefficient  for  basal  slip,  and  since  fi{6  =  0  K)  =  156  GPa  (Table  2)  produces  a 
value  of  energy  per  unit  length  of  basal  plane  edge  dislocations  consistent  with  that  predicted  by  atomic  simulations  (Ta¬ 
ble  8).  For  calculations  of  Section  4.4  of  the  present  paper  that  rely  on  isotropic  shear  and  bulk  moduli  (Table  7),  represen¬ 
tative  values  of  B  =  250  GPa  and  v  =  0.24  have  been  used  (Holm  et  al.,  1999),  leading  to  effective  isotropic  shear  modulus 
C  =  3B(1  -  2v)/(2  +  2v)  =  157  GPa. 
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